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Abstract 

This  report  summarizes  the  FY'92-93  research  activity  at  the  Colorado  Center  for 
Astrodynamics  Research  (CCAR)  on  the  application  of  GPS  based  attitude  determination 
techniques  to  near-Earth  spacecraft.  The  areas  addressed  include  solution  algorithms, 
simulation  of  the  spacecraft  and  environment,  and  ground  testing  of  a  state-of-the-art 
attitude  receiver. 

The  simulation  includes  models  of  a  gravity  gradient  stabilized  gyrostat  in  a  near-Earth, 
J2  perturbed  orbit.  GPS  pseudorange  measurements  for  navigation,  and  phase  difference 
measurements  for  attitude  determination,  were  modeled  including  error  sources  such  as 
Selective  Availability  (SA)  and  multipath.  A  Kalman  Filter  was  designed  to  estimate  the 
position  and  velocity  of  the  vehicle  in  an  Earth  Centered  Inertial  coordinate  frame.  With 
the  simplified  models  used,  positioning  accuracy  was  found  to  be  about  2  meters  RMS 
without  SA  and  30  m  RMS  when  SA  was  included.  For  attitude  determination,  the 
antenna  baselines  were  considered  to  be  1  meter  long  in  an  optimal  orthogonal 
configuration.  The  nominal  performance  of  the  point  (or  instantaneous)  attitude  solution 
was  found  to  be  0.1  deg  RMS.  A  linearized  Kalman  Filter  algorithm  parameterized  by 
Euler  angles  produced  better  performance  at  the  level  of  0  03  deg  RMS.  An  extended 
Kalman  filter  implemented  using  a  quaternion  state  representation  gave  comparable 
performance  when  the  spacecraft  was  within  1  degree  of  local  level,  and  unlike  the  Euler 
angle  version,  was  able  to  perform  equally  well  in  other  spacecraft  orientations.  The 
slowly  varying  multipath  errors  had  a  severe  effect  on  the  filtered  solutions,  degrading  the 
accuracy  to  the  level  of  about  0.1  deg  RMS. 

Ground  experiments  were  performed  with  a  state-of-the-art  GPS  attitude  receiver,  the 
Vector,  manufactured  by  Trimble  Navigation.  This  receiver  has  six  parallel  hardware 
channels,  which  multiplex  among  four  antennas  (1  master  and  3  slaves)  to  produce 
absolute  range  measurements  to  the  master  and  differential  phase  measurements  between 
the  master  and  each  of  the  slaves.  The  antennas  were  mounted  to  one  of  four  test 
structures  to  evaluate  the  measurement  quality  and  the  attitude  determination 
performance.  The  differential  phase  precision  was  found  to  range  from  better  than  2  mm 
1-0  to  as  much  as  1cm  l-o  in  a  particularly  noisy  multipath  environment.  The  Vector 
implements  a  least  squares  algorithm  which  solves  for  vehicle  pitch,  roll,  and  heading 
angles.  The  solution  performance  was  evaluated  in  static  tests  with  the  structure  aligned 
in  a  known  orientation.  The  pointing  accuracy  was  found  to  be  0.2  deg  for  a  1  m  baseline 
structure.  We  also  investigated  an  off-line  self-survey  algorithm  which  establishes  the 
baseline  vectors  in  the  body  fixed  frame  and  determines  the  antenna/receiver  line  biases. 

The  next  step  in  this  project  is  to  test  the  algorithms  on  some  experimental  data  which  we 
recently  received,  collected  from  the  Air  Force  RADCAL  satellite.  In  the  future  we  hope 
to  further  validate  and  improve  our  algorithms  with  other  on-orbit  data  to  improve 
performance,  reliability,  and  robustness  under  varying  conditions.  Our  plans  for  future 
ground  test  are  also  geared  toward  these  goals. 
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1.  Background 

GPS  is  currently  being  described  as  an  “all-purpose”  sensor  for  aircraft  and  spacecraft 
applications.  The  navigation  performance  of  GPS  has  been  widely  demonstrated  on  land, 
marine,  and  airborne  vehicles,  and  to  a  limited  degree  on  near-Earth  spacecraft.  More 
conventional  spacecraft  attitude  sensors  include  conical  Earth  scanners,  magnetometers, 
sun  sensors,  attitude  reference  units  (ARU),  and  star  trackers.  The  highest  performance  is 
obtained  through  a  combination  of  star  tracker  and  ARU,  providing  attitude 
determination  accuracy  better  than  0.01  degrees,  3-a  [Levenson,  et  al.,  1992].  Horizon 
sensors  are  more  coarse,  with  accuracies  at  the  level  of  0.05-0.1  degrees  [Larson  and 
Wertz,  1992].  Instantaneous  GPS  solution  accuracy  is  expected  to  be  in  the  range  of  0.1- 
0.5  degrees,  depending  on  the  vehicle  configuration.  With  dynamic  modeling  and 
measurement  filtering,  further  improvement  is  expected. 

In  the  Iasi  »wo  years,  several  neai  Earth  spacecraft  have  been  launched  which  make  use  of 
GPS  in  some  capacity.  TOPEX/Poseidon  is  a  joint  project  between  the  U.S.  and  France 
to  study  the  dynamics  of  the  Earth's  oceans.  TTie  primary  instruments  are  two  altimeters 
designed  to  measure  the  height  of  the  spacecraft  above  the  ocean  surface.  In  addition, 
TOPEX  carries  a  6  channel  Motorola  GPS  Elemonstration  Receiver  (GPSDR).  This  is  a 
dual  frequency  P-code  and  carrier  tracking  receiver.  Crude  position  solutions  are 
available  on  board  the  satellite;  however,  the  primary  function  of  the  receiver  is  to  collect 
raw  data  for  orbit  post  processing  by  the  Jet  Propulsion  Laboratory  (JPL).  A  recent  paper 
by  Bertiger,  et  al.  [1993]  describes  the  orbit  accuracies  of  3  cm  RMS  radial  and  10  cm 
along  track  and  cross  track  which  are  currently  being  achieved  using  the  TOPEX  data 
together  with  a  network  of  ground  based  GPS  receivers. 

The  Extreme  Ultraviolet  Explorer  (EUVE)  is  also  currently  operating  an  experimental 
receiver  in  a  500  km  altitude  orbit.  The  EUVE  receiver  is  the  engineering  prototype  of 
the  one  to  be  flown  on  TOPEX.  Gold,  et  al.  [1993]  recently  reported  on  orbit 
determination  results  generated  using  GIPSY-OASIS  II  in  both  post- processed  mode 
(with  ground  station  augmentation)  and  in  a  pseudo  real  time  mode.  Unlike  the  TOPEX 
receiver,  the  EUVE  version  only  tracks  LI  C/A  code  and  carrier  signals.  Six  satellites  are 
tracked  at  a  time  from  one  of  two  antennas  mounted  to  the  spacecraft  body.  Their  results 
showed  potential  for  on  board  performance  at  the  several  meter  level  even  in  the  presence 
of  Selective  Availability. 

We  are  aware  of  only  two  tests  of  both  orbit  and  attitude  determination  on  orbit.  The  US 
Air  Force  RADCAL  satellite  is  currently  flying  a  Trimble  Quadrex  receiver  which 
collects  differential  phase  data  from  four  antennas.  Initial  attitude  solutions  have  been 
reported  by  Cohen  and  Lightsey  [1993].  They  have  had  some  difficulty  in  processing  the 
solutions  because  of  unreliable  ground  calibration  of  the  baselines  and  the  severe 
multipath  environment  on  this  spacecraft.  However,  they  are  continuing  to  work  with 
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this  daia.i  An  additional  test  with  the  Trimble  Quadrex  was  recently  performed  on 
Shuttle  Flight  STS-51  [Saunders,  1993],  The  receiver  and  four  antennas  were  carried  in 
the  crew  compartment  and  attached  to  the  orbiter  forward  window  panels.  Results  from 
this  experiment  are  also  expected  within  the  next  several  months. 

Many  upcoming  satellite  missions  have  plans  to  employ  GPS  for  orbit  and  in  some  cases 
attitude  determination,  including  Gravity  Probe  B  [Axelrad  and  Parkinson,  1989], 
Seastar,  Hermes  Columbus  [Munjal,  et  al.,  1992],  the  U.S.  Space  Shuttle  [Madden,  1992], 
and  the  Orbital  Sciences  mobile  communications  satellites,  Orbcomm. 

The  use  of  GPS  for  attitude  determination  on  land,  sea,  and  air  is  somewhat  further 
developed.  There  are  currently  several  commercially  available  attitude  determination 
receivers  (or  systems  of  receivers)  including  those  made  by  Trimble  [Cohen,  1992; 
Trimble  1992],  Texas  Instruments  [Brown,  1990],  Ashtech  [Ferguson,  1991],  and  Adroit 
Systems  [Diefes,  1993].  Many  of  these  were  reported  on  or  demonstrated  at  the  ION 
GPS-93  Meeting  held  in  Salt  Lake  City,  September  22-24, 1993. 

GPS  attitude  determination  systems  have  grown  out  of  the  successes  of  GPS  kinematic 
survey,  wherein  carrier  interferometry  is  used  over  baselines  of  10  km  or  more  to 
precisely  determine  baseline  vectors  from  a  reference  site  to  some  type  of  mobile 
surveyor  on  land,  sea,  or  air.  Kinematic  survey  has  the  potential  to  achieve  accuracies  of 
better  than  10  cm;  however,  it  is  limited  by  multipath  and  range  decorrelation  effects  of 
the  ionosphere,  troposphere,  and  satellite  orbit  errors.  The  fundamental  observation  for 
kinematic  survey  is  a  double  difference  of  carrier  measurements  made  by  the  two 
receivers  to  two  satellites.  The  double  difference  eliminates  both  satellite  specific  errors 
and  receiver  clock  errors.  However,  in  order  to  use  these  observations,  an  integer 
wavelength  ambiguity  must  be  solved  for  each  satellite  pair.  This  generally  requires  extra 
satellites  and  sometimes  may  take  tens  of  minutes  to  complete. 

In  several  ways,  attitude  determination  is  simpler  than  kinematic  survey.  First,  in  some 
attitude  receivers,  a  common  oscillator  is  used  to  make  the  measurements  on  all  antennas, 
thus,  there  is  no  need  to  difference  between  satellites  to  remove  the  receiver  clock  error. 
This  reduces  the  number  of  measurements  required  to  form  a  solution.  Second,  the 
ambiguity  resolution  problem  is  made  easier  because  the  baseline  lengths  and  relative 
orientations  in  the  body  frame  are  fixed  and  known.  The  baselines  are  also  shorter,  thus 
limiting  the  number  of  possible  carrier  cycle  ambiguities  which  must  be  considered. 
Finally,  the  small  separation  between  the  antennas  effectively  eliminates  concerns  about 
ionospheric  and  tropospheric  range  decorrelation.  The  remaining  accuracy  limitations  are 
multipath  reflections,  antenna  phase  center  migration,  antenna  crosstalk,  and  structural 
flexure. 


^  Glenn  Lightsey  bas  provided  us  witb  wbat  be  believes  is  tbe  best  5  hour  data  set  from  RADCAL.  We  will 
begin  work  on  this  witbin  tbe  next  several  weeks. 
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2.  Objectives  and  Benefits 

The  goal  of  our  work  is  to  develop  algorithms  and  investigate  the  potential  for  GPS  based 
attitude  determination  in  near-Earth  spacecraft  applications.  Past  and  current  work 
indicates  that  GPS  has  tremendous  potential  for  providing  real  time  orbit  and  attitude 
information  directly  to  the  vehicle.  This  is  of  interest  for  two  reasons  -  this  data  is 
required  to  both  interpret  data  from  on  board  sensors,  and  as  feedback  to  a  control  system 
to  achieve  the  desired  spacecraft  orbit  and  orientation.  Before  GPS  can  be  used  "in-the- 
loop",  accuracy,  reliability,  and  robustness  of  the  onboard  solutions  must  be  further 
investigated  and  improved. 

The  first  year’s  efforts  on  this  project  were  focused  on  establishing  a  reliable  simulation 
of  the  near-Earth  attitude  determination  performance  and  conducting  ground  experiments 
for  evaluation  of  the  basic  performance  of  GPS.  The  algorithms  derived  in  both  the 
ground  tests  and  simulations  will  soon  be  put  to  an  important  test  in  processing  data 
collected  on  orbit  from  the  RADCAL  mission. 

In  the  simulation  area,  it  was  our  objective  to  set  up  and  validate  a  computer  model  of  the 
spacecraft  environment  including  both  navigation  and  attitude  determination  functions. 
We  would  then  implement  and  compare  various  algorithms,  both  single  point  and  filtered 
solutions,  for  real  time  attitude  determination.  The  simulation  results  have  shown  the 
type  of  performance  that  was  expected  under  nominal  conditions  and  with  measurements 
corrupted  by  Selective  Availability  (for  navigation)  and  multipath.  In  particular,  we 
implemented  a  quaternion  based  estimation  filter  which  is  able  to  perform  under  a  wide 
range  of  spacecraft  configurations. 

In  the  test  area,  we  wanted  to  gain  experience  with  state-of-the-art  GPS  attitude 
determination  hardware  and  real  time  solution  algorithms,  and  to  characterize  both 
measurement  and  solution  quality.  This  is  an  important  step  towards  being  able  to 
process  GPS  measurements  in  real  time  on  orbit 

We  hope  that  this  work  will  provide  important  information  for  designing  attitude 
determination  and  control  systems  for  future  Navy  satellite  missions.  We  expect  to  be 
able  to  suggest  improvements  and  enhancements  of  the  GPS  based  algorithms  which 
should  lead  to  more  accurate,  reliable,  and  robust  on-orbit  performance. 

The  following  section  reviews  the  theory  of  GPS  based  attitude  determination.  Section  4 
describes  the  models  of  the  spacecraft  dynamics,  estimation  algorithms,  and  simulated 
performance.  Section  5  covers  the  algorithms  used  for  ground  testing  and  the 
experimental  results.  Section  6  concludes  the  report  with  a  discussion  of  our  plans  for 
future  work. 


GPS  Based  Spacecraft  Attitude  Determination  9/93 


4 


Colorado  Center  for  Astrodynamcs  Research 


3.  GPS  Attitude  Determination  Overview 

Lisa  M.  Ward 

GPS  based  attitude  determination  is  based  on  observations  of  the  differential  phase 
between  carrier  signals  received  at  two  or  more  antennas.  This  section  discusses  the 
observable  quantities,  how  they  are  obtained  from  GPS  receivers,  and  gives  an  overview 
of  several  attitude  solution  methods. 

3.1  GPS  Signal  and  Observables 

The  GPS  signal  of  interest  for  this  study  is  the  LI  frequency  centered  on  1575.42  MHz. 
Each  GPS  satellite  modulates  this  carrier  with  a  unique  pseudo-random  noise  code  known 
as  the  C/A  code,  at  1.023  MHz  chipping  rate.  The  C/A  code  identifies  the  satellite,  and 
enables  the  receiver  to  measure  the  transit  time  of  the  signal  from  the  satellite  to  the 
receiver.  This  observed  range  based  on  transit  time,  known  as  pseudorange,  is  used  to 
solve  for  the  position  of  the  observer  [c.f.  Milliken  and  Zoller,  1980]. 

The  principal  observable  for  GPS  attitude  determination  is  the  carrier  phase  difference 
between  a  master  antenna  and  one  or  more  slave  antennas.  Figure  3.1  shows  the  incident 
carrier  wave  received  by  a  single  master  and  slave  antenna  pair.  The  difference  in  the 
received  phase,  Atp,  between  the  master  and  the  slave  is  the  actual  measurement 
observable.  The  phase  difference,  A(p,  is  related  to  the  range  difference,  zir,  (both 
expressed  in  cycles)  as  follows, 

Ar  =  A<p  +  k-P+v  (3.1) 

where  k  is  the  integer  number  of  carrier  cycles  in  the  differential  range,  fi  is  the  constant 
fractional  cycle  hardware  delay  between  the  two  antennas,  and  v  is  the  random 
measurement  noise. 


Figure  3.1  GPS  Differential  Phase  Geometry 
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3.2  Relating  Observables  to  Attitude 

As  can  be  seen  in  Figure  3.1,  the  range  difference  is  just  the  projection  of  the  baseline 
vector,  b  (directed  from  the  master  to  the  slave),  onto  the  user-to-GPS  line  of  sight  vector, 
e,  which  gives 


Ar  =  b-e.  (3.2) 

If  we  assume  that  baseline  vector  is  known  in  the  vehicle  or  body  fixed  frame,  and  the 
line  of  sight  vector  is  known  in  an  inertial  (or  Earth  Fixed)  reference  frame,  the 
differential  range  can  be  written  as  a  function  of  the  inertial  to  body  transformation 
matrix  as  follows, 

Ar  =  (b®f®C^e^  (3.3) 

One  representation  of  the  attitude  is  the  transformation  matrix  between  the  I  and  B 
frames. 


3.3  Reference  Frames 

For  this  study  three  reference  frames  are  defined  as  shown  in  Figure  3.2.  The  first  is  the 
inertial  frame,  I.  Here,  is  aligned  with  the  spin  axis  of  the  Earth,  is  in  the  equatorial 
plane  and  points  toward  the  first  point  of  Aries,  and  also  in  the  equatorial  plane  is 
perpendicular  to  i^.  The  next  frame  is  the  orbit  local  frame,  L.  The  vector  points  in  the 
satellite  zenith  or  radial  direction,  is  aligned  with  the  angular  momentum  vector  of  the 
satellite,  and  jL=kL  x  completing  the  orthogonal  set.  Finally,  the  body  frame,  B,  axes 
coincide  with  the  principal  axes  of  the  satellite,  and  are  nominally  aligned  with  the  local 
frame. 

The  transformation  matrix  which  describes  the  rotation  from  the  inertial  to  the  local 
frame  is  derived  from  the  satellite  position.  The  transformation  matrix  describing 
the  rotation  from  the  local  to  the  body  frame  is  derived  from  the  Euler  angles:  yaw,  roll, 
and  pitch.  Here  pitch  is  the  first  rotation  about  the  k®-axis,  roll  is  the  second  rotation 
about  the  j®-axis,  and  yaw  is  the  third  rotation  about  the  i®-axis.  The  total  rotation  from 
the  inertial  to  the  body  frame  is  described  by  the  transformation  matrix 
Instead  of  using  Euler  angles  or  a  transformation  matrix,  quaternions  can  also  be  used  to 
describe  the  rotation  from  between  frames.  The  relationship  between  these  three  attitude 
representations  can  be  found  in  Appendix  A. 
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3.4  Attitude  Solution  Methods 

Deterministic  attitude  determination  is  possible  given  a  minimum  number  of 
observations,  for  example  3  baselines  observing  2  satellites.  However,  there  are  usually 
more  measurements  available  enabling  us  to  optimize  the  solution  in  some  way.  An 
iterative  nonlinear  least  squares  method  is  one  possibility.  Markley  [1988]  describes  a 
fast  least  squares  method  using  singular  value  decomposition  to  find  the  optimal  attitude 
matrix.  Kalman  fillers  are  yet  another  way  to  solve  for  attitude  by  incorporating 
information  about  the  dynamics  into  the  least  squares  process.  An  in  depth  discussion  of 
these  methods  will  be  presented  in  Section  4. 
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4.  Spacecraft  Attitude  and  Orbit  Determination 

Lisa  M.  Ward 

A  simulation  was  developed  to  evaluate  the  performance  of  the  various  attitude  solution 
methods  using  GPS.  The  true  user  satellite  orbit  and  attitude  dynamics  are  modeled  as 
well  as  the  orbit  dynamics  of  the  GPS  satellites.  Measurements  observed  by  the  GPS 
receiver  are  constructed.  The  position  and  velocity  of  the  user  are  estimated  with  the 
navigation  model.  Finally,  the  attitude  of  the  vehicle  is  estimated  using  one  of  the 
solution  methods.  A  description  of  the  simulation  models  and  algorithms  are  presented  in 
Sections  4. 1-4.8.  Section  4.9  provides  specifics  on  the  computer  simulation,  and  4.10 
discusses  the  various  cases  considered  and  the  results  obtained.  The  FORTRAN  code 
implementing  the  simulation  is  available  upon  request. 


4.1  Orbit  Dynamics 

In  this  work,  the  primary  interest  is  in  attitude  determination;  however,  since  the  orbit 
local  frame  as  well  as  the  gravity  gradient  torque  depend  on  the  position  of  the  satellite,  a 
model  of  the  orbital  dynamics  must  be  included  as  well.  A  two-body  orbit  with  oblate 
Earth  perturbations  was  selected  to  describe  the  true  motion  of  the  user  satellite.  These 
equations  of  motion  are. 


where  r  =  +y^  + 

and  where  x,  y,  and  z  are  the  satellite  position  coordinates,  /i  is  the  Earth's  gravitational 
constant,  J2  is  the  Earth's  2nd  zonal  harmonic  coefficient,  and  Re  is  the  radius  of  the 
Earth.  The  position  and  velocity  of  the  user  satellite  are  calculated  by  numerically 
integrating  Equations  4.1. 

The  position  and  velocity  of  the  GPS  satellites  are  calculated  analytically  using  a  simpler 
model  which  includes  only  two-body  effects.  The  orbital  elements  for  the  24  satellite 
"primary"  satellite  constellation  reported  by  Green,  et  al.  [1989]  were  used. 
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4.2  GPS  Orbit  Observations 

Once  the  positions  of  the  user  and  GPS  satellite  are  established,  the  distance,  D,  between 
them  can  be  computed.  TT.i  quantity  that  is  actually  measured  by  a  GPS  receiver  is  the 
pseudorange,  p.  The  p.scudorange  is  equal  to  the  distance  plus  clock  bias  and  noise  errors 
due  to  the  receiver  (bR  and  Up),  bias  errors  due  the  satellite  (.-bs),  delay  through  the 
ionosphere  (e/),  and  Selective  Availability  (s),  as  follows, 

D  —  +  (ygps  ~ y)  ^^gps  ~  ^  ^ 

p  —  D-¥  bff  —  b^  +  £j  +  s  +  Up 


where  Xgps,  ygps,  and  Zgps  are  the  position  coordinates  of  the  GPS  satellite.  The  receiver 
clock  bias  is  not  modeled  in  our  work.  Normally  it  is  solved  for  in  the  navigation  filter 
and  does  not  affect  the  performance.  We  also  ignored  the  bias  errors  due  to  the  satellite 
clock,  orbit,  and  ionosphere. 

Similarly,  the  range  rate,  D,  can  be  computed  as  projection  of  the  relative  velocity  onto 
the  line  of  sight  D  and  the  related  pseudorange  rate,  p  ,  are  given  by. 


^  ~  ~  ~  y^^ygp^  ^^gps  ~  ^^^^gps  ^)]  ^  ^ 

P  =  D  +  //y-/5+i+Up 


.(4. 


where  the  receiver  clock  drift  (fp)  and  satellite  clock  drift  (/jjwill  be  ignored. 

As  a  simple  approximation,  the  measurement  errors  for  both  the  range  ( Up)  and  range 
rate  ( Up)  were  modeled  as  Gaussian  processes  with  zero  mean  and  standard  deviations 
Op  and  Op ,  respectively. 

This  type  of  modeling,  however,  does  not  account  for  Selective  Availability  (SA).  Since 
SA  is  the  largest  error  source  for  non-differential  users  it  is  important  to  include  this 
effect.  A  second  order  ARIMA  model  described  by  Lear  [1992]  is  used  to  more 
accurately  describe  the  range  measurement  noise.  For  each  GPS  satellite  a  selective 
availability  error  is  generated  using 

5*  =  +b,w^+  b^w^_,  (4.4) 

where  Sk  is  the  measurement  error  due  to  selective  availability  at  time  tk  and  wk  is  the 
white  noise  (with  zero  mean  and  unity  standard  deviation)  at  time  tk.  ai  and  bi  are  the 
model  coefficients  for  a  1  second  time  interval  and  are  given  by 
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ai  =  1.993960846  (4.5) 

32  =  -0.9940230677 
bi =  0.998209455 
b2  =  -0.974835034 

Figure  4.2  shows  a  plot  of  the  typical  selective  availability  for  one  GPS  satellite. 


Figure  4.2.  Example  of  Selective  Availability 

The  corresponding  range  rate  measurement  noise.  4.  is  computed  by  approximating  the 
derivative  of  the  selective  availability  with 

(4.6) 

where  At  is  time  interval  between  measurements. 


4.3  Attitude  Kinematics 

The  time  derivative  of  the  inertial  to  body  attitude  matrix  defines  the  inertial  angular 
velocity  vector  of  the  body,  co. 

"C'=(o’‘“C‘  (4.7) 

where  to*'  is  the  cross  product  matrix. 
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0 

CO, 

-CO- 

-CO, 

0 

co, 

CO, 

-co^ 

0 

and  where  (D  =  0),i®+0)2j*  +  cu3k*- 


(4.8) 


Rewriting  this  equation  in  terms  of  Euler  angles  where  63  is  the  first  rotation  about 
3-axis,  02  is  the  second  rotation  about  the  2-axis,  and  61  is  the  third  rotation  about  the 
1-axis  gives 


0,  =Q),-i-(cD2sin0i-»-tU3COS0,)tan02 

02=fi>2cos0,-a)3sin0i  (4.9) 

^3  =(ty2sin^,  -<-ty3Cos^,  )/cos02- 

To  formulate  these  equations  in  terms  of  yaw,  roll,  and  pitch  which  characterize  the  local 
to  body  rotation,  the  orbital  motion  of  the  satellite  must  be  removed.  This  is  done  as 
shown  in  Melvin  [1991]  by  setting  61  equal  to  yaw,  62  equal  to  roll,  and  63  equal  to  pitch 
plus  the  argument  of  latitude,  u. 

In  addition  to  matrices  and  Euler  angles,  quaternions  can  also  be  used  describe  the 
attitude  of  a  body.  Recall  that  the  quaternion  q  is  defined  as 


.^4 


(4.10) 


where 


and 


^4  =co 


a  is  the  rotation  angle  about  the  axis  of  rotation  X.  The  rate  of  change  of  the  quaternion 
is  given  by 


q 


(4.11) 
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where 


^4 

-^3 

92 

^3 

94 

-9i 

-<12 

9\ 

94 

“^1 

-92 

-92 

In  order  for  q  to  represent  the  local  to  body  rotation,  the  O)  vector  in  Equation  4.1 1  must 
be  the  angular  velocity  of  the  body  with  respect  to  the  local  frame,  Since  generally 
the  inertial  angular  velocity  ^tOy  is  needed  for  the  dynamical  equations  below,  the  angular 
velocity  of  the  of  the  local  frame  in  inertial  space  ^co,  must  be  removed  as  follows: 

y,0)f=y,0)yS<  (4.12) 

where  j,a)f="C‘'y^a)f 

and  ^(Oy  =  [0  0  ti]^. 

The  superscripts  denote  the  coordinate  frame  in  which  the  vector  is  expressed. 


4.4  Attitude  Dynamics 

The  general  attitude  motion  of  a  satellite  is  given  by  Euler's  equation: 

— h^  +  ©®xh^  =  N®  (4.13) 

dt 

where  h®  is  the  angular  momentum  vector  of  the  body,  to  is  the  inertial  angular  velocity 
vector  of  the  body,  and  N  is  the  disturbing  torque.  If  the  origin  of  the  body  frame 
coincides  with  the  center  of  mass  of  the  satellite  and  the  axes  are  aligned  with  the 

principal  axes  of  the  satellite  then  h®=[l|0|  l2Q>2  Equation  4.13 

becomes 

I,tu,  +  ty2t03(l3-l3)  =  A(. 

I2t02+ £0,603(1, -l3)  =  (Vj  (4.14) 

^3^3  t0,£02  (I2  ~  Ii )  =  Nj. 

Consider  a  gyrostat  containing  a  momentum  wheel  which  has  its  spin  axis  aligned  with 
the  k-axis  of  the  body  and  has  corstant  angular  momentum  hw  The  total  angular 
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momentum  of  the  satellite  becomes  120)2  130)3  +  and  Equation  4.13 

can  be  rewritten  as  follows, 

I, fi),  + 0)30)3(13 -l2)-h^O)j  =  A^, 

IjtOj+O), 0)3(1, -l3)  +  /i,,0),=iV3  (4.15) 

I3d)3  +  0),0)3(l3-l,)  =N^. 

The  gravity  gradient  torque  Ngg  as  derived  by  [Shuster,  1989]  is  given  by 

NGG=^fx(lf)  (4.16) 

where  r  is  the  satellite  orbit  radius,  f  is  the  unit  vector  in  the  radial  direction,  and  I  is  the 
moment  of  inertia  tensor.  Usually  the  inertia  tensor  is  known  in  the  body  frame.  Since  f 
in  the  local  frame  is  just  [  1  0  0  we  can  easily  convert  it  to  the  body  frame  by  writing 

f  f^=®C^[i  0  of.  (4.17) 

can  be  expressed  in  terms  of  either  quaternions  or  Euler  angles.  If  quaternions  are 
used  the  gravity  gradient  torque  becomes 

4(13  -  h  )(^2^  -  )(^39i  +  9492 ) 

Ngg  =  ^  2(Ii  - 13 ){qi  -^-ql+ql liq^qi  +  q^q2 )  .  (4. 1 8) 

r 

.2(12  ~  II  )(9?  -  92  ~  93  94  )(929l  “  9493 ). 

To  summarize,  the  quaternion  and  inertial  velocity  equations  differential  equations  are 

91  =  [92(^  +  ti)  -  936)2  +  94^] !  2 

92  =[“9l(^^ +“)  +  93^  +94®2l^2 

(4.19) 

93  =  [9lfi>2  “  92^  +  94  (<i^  -  “)1  /  2 

94  =  [~9l^yi  “  92^2  “  93  ~  ")1  /  2 

7,0),  =-(/3-/j)G)2fi)3  +  12£2’(/3-/2)(929,  -9493)(939,  +9^93)  -Ko)^ 

1^6)^  =  -(/,  -  73)0,6)3  +  6n^(7,  -  I^){ql  -  ql  -  q]  +  9^ ){q^q^  +  q^q^ )  +  h^co^ 
I^d)^  =  -(73  -  7,  )6),6)2  +  6fl^(73  -  7,  ){q^  -  ql  -  q]  +  ql  )(q^q,  -  q^q^) 
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where 

If,  on  the  other  hand,  is  written  terms  of  Euler  angles,  substituting  the  gravity- 
gradient  torque  NccCy.r.P)  into  equations  (4.15)  and  neglecting  product  of  angles  and 
eccentricity  gives  the  linearized  Euler  angle  equations  of  motion  derived  by  Melvin 
[1991]. 


=  -[K  +  n(i3  - 12  )]«y  -  K  “  (*i  ^2  - 13)"]^ 

V  =  -[K  +  4n(l3  -  h  )]”'■  +  “  {Ii  +  13)'»]>’  (^*20) 

l3P  =  -3(l2  -I,)«V-l3i‘ 

where  n  is  the  mean  motion  of  the  satellite.  Note  that  the  "r"  in  Equation  4.20  refers  to 
the  roll  angle  and  not  to  the  satellite  radius  as  in  other  equations. 


4.5  GPS  Attitude  Observations 

As  mentioned  earlier,  the  GPS  receivers  measure  a  phase  difference  between  antennas.  If 
the  integer  ambiguity,  k,  and  line  bias,  P,  have  been  resolveo,  '  -'orresponding  range 
difference  can  be  computed  directly  by 

Ar  =  Aq)  +  k  —  p.  (4.21) 

This  range  difference  measurement  Ar  for  a  particular  GPS  satellite  and  baseline  is  the 
projection  of  the  baseline  vector  b  onto  the  user  to  GPS  direction  vector  e  plus  an 
unknown  error,  e,  as  follows: 

Ar  =  b-e+e.  (4.22) 

Since  e  is  measured  in  inertial  coordinates,  it  must  be  transformed  to  the  body  frame  with 
the  inertial  to  body  transformation  matrix.  The  range  difference  is  now  given  as 

Ar  =  {h^f  +  e  (4.23) 

The  error  in  the  range  difference  measurement  for  each  GPS  satellite  is  modeled  as  a 
discrete  first  order  Markov  process  as  follows: 


£k=aek^l+bw,^ 


(4.24) 
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where  £k  is  the  range  error  at  time  4  ,  and  is  the  corresponding  element  of  a  discrete 
white  sequence  with  zero  mean  and  unity  standard  deviation.  The  coefficients  a  and  b  are 
computed  as 


a  =  e-^^,b  =  (l-a^)^a  (4.25) 

where  rris  the  steady  state  standard  deviation  of  £k,  fiis  the  inverse  time  constant,  and  At 
is  the  time  interval.  At  every  time  step  the  error  is  computed  for  each  GPS  satellite  and 
baseline  pair  and  is  added  to  the  range  difference  measurements  for  those  satellites  which 
are  visible.  This  type  of  model  can  be  made  to  represent  slowly  varying  error 
components  such  as  multipath,  line  bias  values,  or  structural  bending.  An  example  of  a 
multipath  type  error  is  illustrated  in  Figure  4.3. 


2000  4000  6000  8000 

time,  sec 


10000 


12000 


Figure  4.3  Multipath  Error  Example 

A  GPS  satellite  is  visible  to  the  user  if  the  cosine  of  the  zenith  angle  is  less  than  the 
sine  of  the  mask  angle  7  or 

cos^  <  cos(90°  -  y)  =  siny  (4.26) 


where  cos^  =  f  •  e . 

See  Figure  4.4  for  an  illustration  of  these  angles.  This  situation  is  ideal  in  that  it  assumes 
that  the  antennas  are  omni-directional;  signals  to  the  antennas  are  only  blocked  by  the 
Earth.  Blockage  due  to  the  satellite  itself  is  not  considered  here.  A  more  accurate  model 
would  specify  a  separate  mask  angle  for  each  anterna  and  take  into  account  satellite 
appendages  which  may  block  the  GPS  signal  to  certain  antennas. 


GPS  Based  Spacecraft  Attitude  Determination  9/93 


15 


Colorado  Center  for  Astrodynamcs  Research 


Figure  4.4.  GPS  Visibility 


4.6  Attitude  Point  Solution  Methods 
Nonlinear  Least  Squares 

Given  m  baselines  and  n  satellites,  an  optimal  attitude  solution  is  obtained  by  minimizing 
the  following  cost  function  with  respect  to  ®C*: 

7(”C')  =  iiwJAr,  -(bff  (*C'  e;)f  (4.27) 

i=l)=l  I-  '  'J 

where  Wy  is  the  measurement  weighting  factor. 

Given  an  initial  estimate  of  the  inertial  to  body  transformation  matrix,  ^Cq,  an  iterative 
approach  may  be  used  to  solve  for  the  best  estimate  of  By  linearizing  the  cost 
function  about  the  current  estimate  of  the  attitude  matrix  a  better  estimate  can  be 
generated  from 


‘C'=  SC  S  (1  +  0") 


(4.28) 


where  1  is  the  identity  matrix  and  0"  is  the  cross  product  matrix  associated  with 
56  =  [56,  S02  5®?]^  as  in  Equation  4.8.  56i,  56  and  S03  represent  small  rotations 
about  the  estimated  body  axes.  Using 

Af,=(b,')"('C;«')  (4.29) 
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to  predict  the  range  difference,  the  cost  function  can  be  rewritten  as 


A'c')=  i iwj8r,j  -(b'f (e”  *c; e')] 

(4.30) 

where 

Sr  =  At  —  Af  . . 

>;  V  v 

(4.31) 

The  second  term  in  the  cost  function  can  be  rearranged  as  follows 

(bf  f  (©*  ej)  =  ejf  BrSe  (4.32) 


where  B)'  is  the  cross  product  matrix  associated  with  baseline  vector  bi. 
Then  the  linearized  cost  function  can  be  written  as 


^(se)| 


*ci 


where  the  measurement  gradient,  H,  is  formed  by 


(4.33) 


(4.34) 


Now,  to  find  the  best  estimate  of  solve  the  system  of  equations  H,-,50  =  5r,y  for  60. 

Recompute  the  attitude  matrix  using  Equation  4.28  and  repeat  until  50  gets  sufficiently 
small. 

The  quaternion  formulation  is  similar  except  the  linearized  version  of  the  attitude  matrix 
in  Equation  4.28  takes  on  the  form 

={l  +  2Q^fc'o  (4.35) 


where  is  the  cross  product  matrix  associated  with  dq  =  [5g,  8^2 
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SVD  Point  Solution 

Markley  [1988]  derives  a  fast  method  for  solving  a  similar  cost  function  using  singular 
value  decomposition.  This  cost  function,  known  as  Wahba’s  problem,  is  given  as 


i=l 


(4.36) 


where  the  unit  vectors  sf  and  s/  are  the  direction  vectors  to  the  /th  object  in  the  body 
frame  and  inertial  frame,  respectively.  Cohen  [1992]  shows  how  the  Markley  solution 
can  be  applied  GPS  attitude  determination.  One  advantage  of  this  method  is  that  it  does 
not  require  any  prior  knowledge  of  the  attitude.  Hence,  this  solution  could  be  used  in 
initializing  a  Kalman  filter. 

To  derive  the  SVD  solution  (following  Cohen  [1992])  convert  the  cost  function  in 
Equation  4.27  to  matrix  form 

7(C)  =  ||wi,'^(AR  -  C  (4.37) 


where  AR  is  the  differential  range  measurements  assembled  into  a  (mxn)  as  follows 


AR  = 


Ar,i  Ar,2 

Ar2i  Ar22 


Ar, 


In 


Ar, 


ml 


Ar, 


mw  J 


(4.38) 


and  where  C  is  the  estimate  of  the  inertial  to  body  transformation  matrix,  B  are  the  m 
baseline  vectors  in  the  body  frame  concatenated  into  a  (3xm)  matrix,  E  are  the  n  user  to 
GPS  direction  vectors  in  the  inertial  frame  concatenated  into  a  (3xn)  matrix,  and  Wb  and 
We  are  weighting  matrices  for  the  baseline  and  direction  vectors,  respectively.  The 
superscripts  have  been  omitted  here  for  simplicity. 

The  Frobenius  norm  for  a  matrix,  M,  is  defined  as 

||Mlg  s  trace(M‘^  m).  (4.39) 

Using  the  trace  equivalence  of  the  Frobenius  norm,  the  matrix  form  of  the  cost  function 
can  be  expanded  to  give 
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y(c)  =  trace^W^  Wg  AR  j 

-  2  trace^^W^  ^r'*'  Wg  C  E  j  (4.40) 

+  trace^^W^  (c  e)^  B  Wb  C  E 

Because  the  first  term  is  independent  of  attitude,  it  drops  out  when  minimizing  with 
respect  to  the  attitude.  With  the  proper  choice  of  Wb  the  last  term  also  becomes 
independent  of  attitude  and  drops  out  Choosing  Wb  will  be  discussed  below. 

After  rearranging  the  second  term,  the  trace  being  invariant  to  cyclic  permutations,  the 
problem  now  is  one  of  maximizing  the  new  cost  function 

y'(c)  =  trace(c  E  AR"^  Wg  B"^)  s  trace(c  G"^)  (4.41 ) 

where  G  5  B  Wg  AR  Wg  E^. 

The  f'ptimal  value  of  the  attitude  matrix  can  be  found  by  taking  the  singular  value 
decomposition  of  the  matrix  G.  The  singular  value  decomposition  is  given  by 

G  =  U  S  (4.42) 

where  U  and  V  are  orthogonal  matrices  and  S  is  a  diagonal  matrix  with 

•^11  -  ^22  -  *^33  -  0- 

To  ensure  that  the  optimal  attitude  matrix  is  right-handed,  the  signs  of  the  third  column  of 
U  and  V  matrices  are  adjusted  so  that  their  determinant  is  equal  to  one.  These  new 
matrices  will  be  called  U+  and  \+.  The  sign  of  the  last  element  in  the  S  matrix  is  also 
adjusted  so  that  the  new  matrix.  S',  will  have 

5^'3=det(U)det(V)533.  (4.43) 

Now  G  can  be  rewritten  as 


G  =  S'  \l. 


(4.44) 


and  the  cost  function  becomes 
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nC)  =  trace(c  S'  Uj)  =  trace(s'  Uj  C  V^) 
=  trace(S'A) 


(4.45) 


where  AsUjCV^. 


(4.46) 


Since  A  is  an  orthogonal  matrix,  an  Euler  angle/axis  representation  can  be  used  to  express 
it  as  follows, 

cosa  +  X^(l-cosa)  XjSina  +  XjXjCl-cosa)  -XjSina  +  Cl-cosalXjXj 

A=  -X,3sina  +  X.2Xj(l-cosa)  cosa  +  X^fl-cosa)  A,,sina  +  (1  -  cosalX^Xj  (4.47) 

Xjsina  +  XjXjfl-cosa)  -XiSina  +  XjXjCl-cosa)  cosa  +  (l-cosa)A,3 

where  X  =  [X.,  X2  X3]^  is  the  rotation  axis  and  a  is  the  rotation  angle.  Substituting 
this  into  Equation  4.45  gives 

y'(C)  =  5,'[cosa  +  X^j(l  -  cosa)] 

+  S^[cosa  +  X2(l-cosa)]  (4.48) 

+  S3[cosa  +  X3(l  -  cosa)]. 

Since  5,',  S2,  53  >  0  and  X-^,  X-^  ^  1,  then  7'(C)  is  maximized  for  a  =  0.  A  then 

A 

is  equal  to  the  identity  matrix,  and  the  optimal  value  of  C  is  given  by 

U.Vl.  (4.49) 


Cohen  [1992]  points  out  that  to  realize  the  full  potential  of  the  SVD  solution,  the 
baselines  should  form  an  orthogonal  set.  Wb  is  a  measure  of  how  much  the 
measurements  have  to  be  scaled  to  transform  the  baselines  into  an  equivalent  orthonormal 
basis.  So  Wb  is  chosen  such  that 

Wb  =  VbSb'v5  (4.50) 

where  the  singular  value  decomposition  of  B  is  given  by 

B  =  UbSbV5.  (4.51) 

and  S’2  denotes  taking  the  reciprocal  squared  of  each  diagonal  element.  Using  this  SVD 
method  a  new  point  solution  may  be  formed  at  each  measurement  epoch. 
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4.7  Attitude  Filters 

Incorporating  known  dynamics  into  the  attitude  determination  process  will  allow  us  to 
obtain  a  better  solution  than  using  the  solutions  presented  above.  As  [Leffens,  1982] 
points  out,  the  Kalman  filter  does  this  by  using  a  dynamical  model  to  propagate  the 
attitude  in  time  along  with  a  measurement  model  to  obtain  the  best  estimate  of  the 
attitude.  For  this  study  two  versions  of  a  Kalman  filter  have  been  implemented  to 
estimate  spacecraft  attitude.  One  is  an  Euler  angle  formulation,  while  the  other  is  a 
quaternion  formulation.  The  general  equations  for  a  Kalman  filter  are  presented  next 
followed  by  the  specific  equations  for  each  filter. 

General  Kalman  Filter  Equations 

(For  background  information  see,  for  example  Gelb  [  1984]). 

The  time  update  of  the  Kalman  filter  is  given  by 


**  =  4-1+ 

/*-i  (4.52) 

where  Xk  and  Pk  are  the  state  vector  and  covariance  matrix  at  time  tk-  The  superscripts,  - 
and  +,  denote  prior  and  post  measurement  update,  respectively  .  d>  is  the  state  u-ansition 
matrix  from  time  tk-i  to  time  4^,  and  Q  is  the  process  noise  matrix. 

The  state  transition  matrix  satisfies  the  equation 


<i>  =  Fd>,  <I)(ro)  =  l, 


(4.53) 


where  F  =  dx  /  9x  is  the  linearized  dynamics  matrix.  An  approximate  solution  for  <1)  is 
given  by  Melvin  [1991].  The  result  is 

<b,j=l^  +  &,{Fy+&l(i,F^  +  St(i^FD))  (4.54) 

where  is  the  identity  matrix  and  At  is  the  filter  time  increment 

The  measurement  update  equations  for  the  standard  Kalman  Filter  are 
K,=P;HJ(H.P.-H;  +  R.r 

X*  =  xj  +  Kj  Az*  (4.55) 

P;  =  (l-K*  HjP; 
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where  K  is  the  Kalman  gain,  R  is  the  variance  of  the  measurement  noise,  and  Az  is  the 
measurement  residual. 

Euler  Angle  Filter 

The  state  vector  for  the  Euler  angle  filter  contains  the  yaw,  roll,  pitch,  and  their  rates, 
(y,  r,  p,  y,  r,  p) .  The  linearized  equations  of  motion  (4.20)  are  numerically  integrated  to 
advance  the  state  estimate  to  the  current  time.  The  covariance  matrix  is  propagated 
forward  with  the  state  transition  matrix.  The  dynamics  matrix  used  to  compute  the  state 
transition  matrix  is  derived  directly  from  Equations  4.20  and  is  shown  in  Appendix  B. 

For  each  measurement  taken  at  the  current  time  the  user  to  GPS  direction  vectors  are 
converted  from  the  inertial  to  the  local  frame  (using  the  position  estimated  by  the 
navigation  filter).  Since  the  Euler  angles  cannot  be  simply  added  together,  the  filter  must 
always  linearize  about  the  local  level.  Thus,  the  residual  (5ry  -  for  baseline  /  and  satellite 
j,)  used  for  the  measurement  update  is  calculated  according  to  Equation  (4.31)  with  the 
Ar.j  given  by 


Af^=(b,Y(l  +  0")e,^  (4.56) 

The  501,  502.  and  503  used  to  form  0*  are  the  yaw,  roll,  and  pitch  angles,  respectively. 
The  measurement  gradient  vector  for  the  angle  states,  H,y,  for  a  single  baseline  /  and 
satellite  j,  then  is 

H,^=(e{-)^Br.  (4.57) 

and  the  elements  corresponding  to  the  rate  states  are  zero.  The  measurement  update  takes 
the  form  shown  in  Equation  4.55. 

One  disadvantage  to  this  formulation  is  that  deviations  from  the  local  level  must  remain 
small  for  the  filter  to  converge  since  we  are  linearizing  about  the  local  frame.  Moving  to 
a  quaternion  filter,  however,  alleviates  this  problem. 

Quaternion  Filter 

The  state  vector  for  the  quaternion  filter  contains  q\,q2,q'i,  and  the  inertial  angular 
velocity  of  the  body.  The  fourth  element  of  q,  <74,  is  included  in  the  state  propagation,  but 
is  not  used  in  the  measurement  update  or  covariance.  This  is  based  on  Lefferts,  et  al. 
[1982]  and  is  described  further  in  the  following. 

Similar  to  the  Euler  angle  filter,  the  quaternion  equations  of  motion  (4.19)  are 
numerically  integrated  to  advance  the  state  estimate  to  the  current  time,  and  the 
covariance  matrix  is  propagated  forward  with  the  state  transition  matrix.  The  elements  of 
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the  dynamics  matrix  used  to  compute  the  state  transition  matrix  are  presented  in 
Appendix  B. 

I 

The  best  estimate  of  the  transformation  matrix  from  inertial  to  body  coordinates  Cq  is 
computed  based  on  the  quaternion  state.  The  range  difference  prediction  is  then, 

Afy=(bf)'^  ®ciej.  (4.58) 

and  the  measurement  residual  is  computed  by  subtracting  this  from  the  observed  Ary 

If  we  linearize  about  the  current  attitude  estimate,  an  improved  estimate  of  the 
transformation  matrix  can  be  modeled  as, 

®C^=(i  +  2Q'')®c4  (4.59) 

where  is  the  cross  product  matrix  for  a  small  correction  rotation  represented  by  the 
quaternion  5q.  Thus  the  measurement  gradient  matrix,  H,  for  the  quaternion  states  is 

Hy  =2(*C^eJ)'^  Bf.  (4.60) 

The  elements  of  the  gradient  matrix  corresponding  to  the  velocity  states  are  zero. 

The  Kalman  gain  calculation  and  covariance  update  follow  Equations  4.55,  but  the  state 
update  is  different.  The  first  three  elements  of  the  correction  quaternion  5q  are  computed 
by  the  filter.  The  fourth  element,  5q4,  is  set  to  normalize  the  quaternion  as  follows, 

5^4  .  (4.61) 

The  correction  quaternion,  6q,  is  then  combined  with  the  apriori  estimate,  according 
to  the  quaternion  composition  rule.  The  general  form  of  the  composition  is  given  by  [c.f. 
Shuster], 


q"  =  q'®q 


'  <l4 

93 

-92 

9i'' 

'9r 

-^3 

94 

9i' 

92 

92 

92 

-9i' 

94 

93 
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.-9i' 

-92 

-93 

94. 

.94. 

(4.62) 


where  q  is  the  first  rotation  expressed  in  the  coordinates  of  the  original  frame;  q'  is  the 
second  rotation  expressed  in  the  coordinates  of  the  intermediate  frame;  and  q"  is  the 
combined  rotation  expressed  in  the  final  frame. 

In  our  case,  q]^ ,  is  the  first  rotation,  and  5q  is  the  second. 
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4.8  Navigation  Filter 

The  attitude  filter  needs  information  about  the  user  satellite's  position  to  calculate  the 
inertial  to  local  transformation.  The  gravity  gradient  torque  also  depends  on  the  user's 
position.  Usually  only  an  estimate  of  the  position  is  available,  so  an  extended  Kalman 
filter,  like  that  used  in  the  attitude  filter,  was  considered.  The  state,  containing  the 
position  and  velocity,  is  propagated  by  integrating  the  orbit  equations  of  motion  (4.1) 
while  the  covariance  matrix  is  propagated  with  the  state  transition  matrix.  Appendix  B 
contains  the  non-zero  elements  of  the  dynamics  matrix  used  to  compute  the  state 
transition  matrix. 

For  each  pseudorange  measurement,  p,  the  residual,  5p,  is  given  by 

^P=p-\^Kl\  (4.63) 

where  the  current  estimate  of  the  relative  distance  between  the  GPS  and  user  satellite  is 

-^0 

%  (4.64) 

A 

20_ 

and  where  Xq,  and  Zq  are  the  current  estimates  of  the  position  coordinates.  The 
measurement  gradient  matrix  can  be  computed  as  follows, 

H=[e,  62  ej  0  0  O]  (4.65) 

where  e  =  4^. 

Similarly,  the  range  rate  measurement  residual,  6p  ,  is  given  by 

5p  =  p-v„,e  (4.66) 


where  the  current  estimate  of  the  relative  velocity  between  the  GPS  satellite  and  user  is 


(4.67) 


and  where  Xq,  yQ,  and  Zq  are  the  current  estimates  of  the  velocity  coordinates.  The 
measurement  gradient  vector  can  be  computed  as 
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4.9  Simulation 

To  study  various  attitude  determination  methods  a  simulation  was  developed 
implementing  the  models  described  in  the  previous  sections.  It  is  comprised  of  three 
main  parts  --  the  truth  models  including  the  user  attitude  and  orbital  dynamics,  the 
measurement  models  for  observed  range,  range  rate,  and  range  differences,  and  finally 
the  estimation  models  for  both  attitude  and  navigation.  The  main  modules  and  data  flow 
within  the  simulation  are  illustrated  in  Figure  4.5. 
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Care  was  taken  to  separate  each  model  or  major  calculation  into  different  modules. 
Coding  in  this  manner  affords  a  number  of  advantages.  For  example,  changes  and 
upgrades  can  be  implemented  with  ease  by  just  replacing  the  old  subroutine  with  the  new 
one.  In  addition,  each  subroutine  can  be  independently  tested  which  simplifies 
debugging  and  verification. 

Initial  orbit  and  attitude  conditions,  satellite  configuration  data,  and  filter  parameters  are 
defined  in  an  input  file.  These  values  can  then  be  easily  be  changed  to  see  the  effect  that 
various  parameters  have  on  the  results.  For  this  study  a  nominal  set  of  input  parameters 
for  the  nadir  pointing  satellite  are  shown  in  Table  4.1. 
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Table  4. 1  Nominal  Simulation  Initialization  Parameters 
Orbital  elements  for  the  user  satellite 


semimajor  axis:  7000  km 

eccentricity:  0.0 

inclination:  60  deg 

ascending  node:  135  deg 

argument  of  perigee:  292  deg 

mean  anomaly:  345  deg 


Satellite  Configuration 

Moments  of  inertia:  Ii  =  1000  kg  m^ 

12  =  2000  kg  m2 

13  =  2010  kg  m2 

Angular  momentum  wheel:  hw  =  SO  kg  m2/s 

Baselines:  3  baselines  along  principal  axes  of  body,  each  1  meter  long 


Measurement  Errors 

Range  difference:  0  =  2  mm 

inverse  time  constant  P=  100/sec  (-white  noise) 
Pseudorange:  o  =  6  m  (white  noise) 

Pseudorange  rate:  o  =  ,025  m/s  (white  noise) 

Initial  attitude 

yo  =  rQ  =  pQ=0.5  deg 

True: 

yo  =  tq  =  po  =  0-001  deg  /  s 


>'0  =  ''0  =  Po=Odeg 

Estimated: 

:yo  =  ^o  =  Po=o  deg/s 

Initial  position 

True:  (calculated  from  orbital  elements  above) 

Estimated:  (expressed  as  a  deviation  from  the  truth) 

Sxo  =  5yo  =  5zo=lm 
Sxq  =  Sy^  =  Szo=  0.001  m  /  s 


Quaternion  filter  parameters 

quaternion  state  process  noise  =lxl0'9  (rad)2 
inertial  velocity  state  process  noise  =1x10*12  (rad/s)2 
measurement  noise  variance  =  4x10*^  m2 

Euler  angle  filter  parameters 

angle  state  process  noise  =  IxlO*^  (deg)2 
angle  rate  state  process  noise  =  1x10*^  (deg/s)2 
measurement  noise  variance  =  4x10*^  m2 
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Each  simulation  was  run  for  200  minutes  or  approximately  2  orbital  periods.  The  time 
between  measurement  updates  was  10  seconds.  A  variable  step  integrator  was  used  to 
propagate  the  state  between  epochs.  This  integrator  is  a  4th  order  Runge-Kutta  library 
routine. 

The  SVD  point  solution,  the  Euler  angle  filter,  and  the  quaternion  filter  were  all  executed 
with  the  above  parameters.  Modifications  were  then  made  to  this  nadir  pointing  scenario 
to  illustrate  the  performance  of  various  algorithms.  For  example,  the  initial  attitude  was 
changed  so  that  the  user  satellite  was  not  nadir  pointing.  Again,  all  three  attitude 
determination  methods  were  executed.  Other  cases  studied  include  adding  a  multipath 
type  error  to  the  range  difference  measurements,  defining  a  larger  mask  angle,  and  using 
a  baseline  configuration  that  is  not  optimal.  In  addition,  the  results  from  the  navigation 
filter  are  presented  below  using  1)  random  measurement  noise,  and  2)  selective 
availability. 


4.10  Results  and  Discussion 

Summaries  of  the  attitude  determination  and  navigation  results  are  presented  below  in 
Tables  4.2  and  4.3  respectively.  Cases  using  the  SVD  point  solution  are  denoted  with 
SVD,  those  using  the  Euler  angle  filter  are  denoted  with  E,  and  those  using  the  quaternion 
filter  are  denoted  with  Q. 


Table  4.2  Attitude  Detetmination  Results 


Case 

Method 

RMS  of  Errors  in  degrees  x  10*2 

yaw 

roll 

pitch 

3-axis 

Nadir  Pointing 

SVD 

4.76 

4.85 

4.78 

8.31 

Nadir  Pointing 

E 

1.35 

1.04 

1.89 

2.55 

Nadir  Pointing 

Q 

1.26 

1.16 

1.62 

2.36 

Off-Nadir  Pointing 

SVD 

4.78 

4.85 

4.79 

8.32 

Off-Nadir  Pointing 

E 

6.60 

6.09 

6.12 

10.63 

Off-Nadir  Pointing 

Q 

1.25 

1.16 

1.62 

2.35 

Multipath  Noise 

Q 

4.28 

5.62 

5.19 

8.76 

Multipath  Noise 

SVD 

4.35 

5.56 

5.15 

8.77 

Non-optimal  Baselines 

Q 

1.85 

1.25 

2.29 

3.20 

Equatorial  Orbit 

Q 

1.25 

1.06 

1.50 

2.22 

Polar  Orbit 

Q 

1.22 

1.13 

1.67 

2.36 

Larger  Mask  Angle 

Q 

2.11 

1.75 

2.39 

3.63 

Selective  Availability 

_Q _ 

1.10 

1.08 

1.11 

1.90 
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Table  4.3  Navigation  Filter  Results 


Case 

RMS  of  Position  Errors  in  meters 

X 

y 

z 

3-axis 

Random  Noise 

0.70 

0.69 

0.66 

1.19 

BHBBWWIW 

17.84 

17.01 

17.91 

30.47 

Case 

RMS  of  Velocity  Errors  in  m/s 

X 

y 

z 

3-axis 

Random  Noise 

7.35x10-3 

5.14x10-3 

3.89x10-3 

9.78x10-3 

Nadir  Pointing  Case 

The  true  yaw,  roll,  and  pitch  generated  by  the  nadir  pointing  case  described  above  are 
shown  in  Figure  4.6.  A  plot  of  the  SVD  point  solution  errors  are  presented  in  Figure  4.7. 
Since  errors  for  all  three  of  the  Euler  angles  are  similar,  the  root-sum-square  (RSS)  of  the 
errors  has  been  computed  to  give  a  3-axis  attitude  error.  The  root-mean-square  (RMS)  of 
the  3-axis  errors  for  the  SVD  point  solution  is  0.084  degrees. 
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Figure  4.7  Nadir  Pointing  SVD  Point  Solution  Errors.  RMS  =  0.084  deg 

With  a  Kalman  filter,  however,  this  error  can  be  improved  since  knowledge  about  the 
dynamics  has  been  included.  As  we  see  in  Figure  4.8,  implementing  the  Euler  angle  filter 
lowers  the  RMS  by  a  factor  of  three  to  0.026  degrees.  The  quaternion  filter  results  for 
this  case  are  shown  in  Figure  4.9.  (Although  quaternions  are  used  in  the  filter 
calculations,  they  have  been  converted  to  Euler  angles  on  output  for  easier  comparison.) 
As  expected,  the  results  are  very  similar  to  those  generated  by  the  Euler  angle  filler,  the 
RMS  being  0.0236  degrees. 
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Figure  4.8  Nadir  Pointing  Euler  Angle  Filter  Errors.  RMS  =  0.025  deg 
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Figure  4.9  Nadir  Pointing  Quaternion  Filter  Errors.  RMS  =  0.0236  deg 


OfT-Nadiir  Pointing  Case 

The  quaternion  formulation  is  better  than  the  Euler  angle  formulation  in  an  off-nadir 
pointing  scenario.  For  example.  Figure  4.10  shows  the  true  attitude  motion  using  initial 
yaw,  roll,  and  pitch  angles  of  3  degrees.  With  this  scenario  the  RMS  of  the  Euler  angle 
filter  3-axis  errors  is  0.1064  degrees,  while  the  RM  S  of  the  quaternion  filter  remains  at 
0.0235  degrees.  The  3-axis  errors  resulting  from  this  off-nadir  pointing  case  are  shown  in 
Figures  4. 1 1  and  4. 12. 
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Figure  4.10  True  Attitude  Off-Nadir  Pointing  Case 
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Figure  4. 11  Off-Nadir  Pointing  Euler  Angle  Filter  Errors.  RMS  =  0.1064  deg 
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Figure  4. 12  Off-Nadir  Pointing  Quaternion  Filter  Errors.  RMS  =0.0235  deg 
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Multipath  Case 

In  the  next  case  a  multipath  type  error  was  added  to  the  range  difference  measurements 
instead  of  adding  only  random  noise.  The  inverse  time  constant  was  adjusted  to  1/300 
seconds**  corresponding  to  a  time  constant  of  5  minutes.  The  RMS  of  the  errors 
increases  to  0.0876  degrees.  Furthermore,  in  Figure  4.13  it  is  evident  that  the  error  is 
more  structured  in  the  multipath  case  than  in  the  random  noise  case.  In  fact,  the  filter 
does  not  work  much  better  than  the  SVD  point  solution  as  shown  in  Figure  4.14.  This 
indicates  that  using  this  filter  in  a  multipath  environment  has  diminishing  returns,  and  a 
better  way  of  dealing  with  multipath  is  needed. 
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Figure  4. 13  Multipath  Environment  Quaternion  Filter  Errors.  RMS  =  0.0876  deg 
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Figure  4.14  Multipath  Environment  SVD  Point  Solution  Errors.  RMS  =  0.0877  deg 
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Non-Optimal  Baselines  Case 

Another  factor  which  can  degrade  the  performance  of  the  attitude  filter  is  poor  antenna 
placement.  In  the  nadir  pointing  case  above  the  baselines  are  optimally  placed  in  an 
orthogonal  configuration.  If  the  baselines  are  not  optimally  configured,  the  pertormance 
of  the  filter  changes.  Consider  the  baseline  vectors  illustrated  in  Figure  4.15  where 
bi  =  [1  0  0]T,  b2  =  [0  0  -0.5]’^,  and  b3  =  [0  0  0.5F.  Using  this  baseline  configuration  we 
would  expect  the  roll  performance  to  be  the  best  because  it  is  observed  by  two  baselines  a 
1  m  and  an  0.5  m  baseline.  The  next  best  should  be  yaw  which  is  observed  by  the  two 
0.5  m  baselines.  Finally  the  pitch  axis  is  only  seen  by  the  long  1  m  baseline.  Figure  4.16 
confirms  this  result.  The  RMS  of  the  roll  errors  0.0125  deg,  while  tlie  RMS  for  the  yaw 
and  pitch  errors  are  0.0185  and  0.0229  ,  respectively.  This  was  not  a  particularly  bad 
configuration,  but  other  more  severe  cases  are  certainly  possible. 


Figure  4. 15  Non-Optimal  Baselines 
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Figure  4.16  Non-optimal  Baseline  Configuration 
RMS(y)=  0.0185  deg;  RMS(r)  =  0.0125  deg;  RMS(p)  =  0.0229  deg 
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Larger  Mask  Angle  Case 

The  number  of  visible  GPS  satellites  also  affects  the  accuracy  of  the  attitude  solution. 
For  instance,  if  the  mask  angle  is  increased  to  15  degrees  blocking  off  more  of  the  GPS 
satellites,  the  RMS  of  the  3-axis  error  increases  to  0.0363  degrees.  We  can  see  in  Figure 
4. 17  that  although  spread  of  errors  is  greater  than  when  the  mask  angle  is  zero,  the 
character  of  result  does  not  change. 


0.10 

0.08 


time,  sec 


Figure  4, 17  Larger  Mask  Angle.  RMS  =  0.0363  deg 

Navigation  Filter  with  Random  Noise 

The  results  from  the  navigation  filter  for  the  nadir  pointing  case  using  random  range 
measurement  errors  are  presented  in  Figures  4.18  and  4.19.  The  x,  y,  and  .  Tors  are 
have  been  combined  to  give  a  3-axis  RSS  position  error  and  the  x,  y,  and  z  errors  have 
been  combined  to  give  a  3-axis  velocity  error.  In  this  particular  case  the  velocity 
measurements  are  better  than  the  position  measurements.  Since  the  velocity  and  position 
states  are  coupled,  the  filter  can  rely  more  on  the  velocity  measurements  resulting  in  the 
apparent  coloring  of  the  position  emors.  The  RMS  of  the  position  and  velocity  emors  are 
1.19  meters  and  9.78  x  10*^  meters  per  second,  respectively. 
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Figure  4.18  Random  Noise  Navigation  Filter  Position  Errors.  RMS  =  1.19  m 
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Figure  4.19  Random  Noise  Navigation  Filter  Velocity  Errors.  RMS  =  9.78  x  10'^  m/s 


Navigation  Filter  with  Selective  Availability 

If  selective  availability  is  included  in  the  range  measurement  noise,  the  performance  of 
the  navigation  filter  worsens  considerably.  The  average  standard  deviation  of  the  errors 
for  this  case  jumps  to  30.5  meters,  as  demonstrated  in  Figure  4.20.  But,  even  though  the 


GPS  Based  Spacecraft  Attitude  Determination  9/93 


39 


Colorado  Center  for  Astrodynamics  Research 


position  error  is  relatively  high,  the  attitude  solution  can  still  be  determined  with  accuracy 
as  indicated  by  an  RMS  of  the  attitude  errors  of  0.019  degrees.  To  accommodate  the  1 
second  interval  ARIMA  model,  the  measurement  update  interval  here  was  changed  to  1 
second. 


Figure  4.20  Navigation  Filter  Selective  Availability  Errors.  RMS  =  30.47  m 


Conclusions 

Nominal  instantaneous  attitude  accuracy  was  found  to  be  at  the  level  of  0.08  deg  RMS 
for  all  three  angles,  based  on  the  measurement  models  which  we  used.  A  factor  of  3 
improvement  was  gained  by  using  a  Kalman  Filter  in  either  implementation.  The 
quaternion  form  was  able  to  maintain  this  accuracy  in  off-nominal  attitudes;  however, 
thus  far  we  have  not  considered  changes  in  the  satellite  visibility  due  to  changing  vehicle 
attitude.  This  is  certainly  an  area  for  future  improvements. 

The  attitude  filter  solution  degradation  in  the  presence  of  even  small  multipath  type 
reflections  is  worse  than  expected.  In  fact,  it  appears  that  there  is  no  benefit  to  filtering  in 
this  case.  We  will  continue  to  investigate  both  the  assumptions  and  the  models  used,  to 
make  sure  that  the  problem  is  accurately  represented,  and,  if  so,  to  try  to  develop 
techniques  to  improve  performance. 
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5.  Algorithms  and  Ground  Testing 

Bruce  C.  Chesley 

This  section  describes  algorithms,  testing,  and  results  for  GPS  based  attitude 
determination  on  the  ground.  The  algorithms  and  experiments  address  static  attitude 
solutions,  baseline  survey,  and  integer  ambiguity  resolution.  Before  describing  the 
ground  algoritms  and  experiments,  a  summary  of  GPS  attitude  receivers  and  error  sources 
is  given.  Note  that  many  of  the  test  descriptions  and  results  have  been  previously 
presented  by  Axelrad  and  Chesley  [1993]. 


5.1  Attitude  Receiver  Comparison 

There  is  currently  a  great  deal  of  interest  in  GPS  based  attitude  determination  in  the 
navigation  community.  Several  experiments  have  been  conducted  utilizing  somewhat 
different  approaches.  The  hardware  configurations  can  be  divided  into  two  catagories  - 
those  using  multiple  independent  receivers  (usually  driven  from  a  common  oscillator), 
and  implementations  of  a  single  receiver  with  multiple  antennas.  The  single  receiver 
versions  have  varying  degrees  of  hardware  integration. 

Independent  receivers  -  Currently,  the  most  widely  reported  approach  to  GPS  based 
attitude  determination  is  to  use  two  or  more  separate  receivers  to  collect  data  from 
multiple  antennas  [c.f.  Martin-Niera  and  Lucas,  1992;  Lu,  et  al.,  1993;  Brown  and  Evans, 
1990;  Diefes,  1993].  The  phase  measurements  from  the  receivers  are  processed  together 
to  eliminate  common  mode  errors  and  solve  for  the  antenna  array  orientation.  Generally 
this  requires  an  additional  satellite  to  resolve  the  oscillator  offsets  between  the  receivers 
or  to  permit  measurement  double  differencing  to  eliminate  its  effect. 

The  use  of  a  common  reference  oscillator  for  all  receivers  may  eliminate  the  need  for  an 
additional  observation.  However,  depending  on  the  design,  the  relative  offset  between 
the  receivers  is  still  subject  to  hardware  dependent  delays  which  can  be  a  function  of  the 
receiver  power  up  state  or  may  vary  due  to  differences  in  temperature  or  supply  voltage. 
Thus,  even  in  a  common  oscillator  configuration,  double  differences  are  generally  used 
for  multiple  receiver  attitude  determination. 

Ash  tech  -  The  approach  taken  in  the  Ashtech  3DF  receiver  is  to  devote  a  separate 
hardware  channel  for  each  antenna  baseline  driven  with  a  single  oscillator.  The  hardware 
has  been  integrated  into  a  single  package,  but  the  signal  processing  approach  is 
essentially  the  same  as  the  independent  receiver  approach  described  above.  Sample 
testing  of  this  receiver  has  been  documented  by  Ferguson  [1991],  and  Van  Graas  and 
Braasch  [1991]. 

Texas  Instruments  -  TI  conducted  several  development  tests  of  an  attitude  determination 
receiver  that  were  reported  by  Brown  and  Ward  [1990]  and  Brown  and  Evans  [1992]. 
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These  tests  used  a  modified  AN/PSN-9  receiver,  an  off-the-shelf  production  model  from 
TI.  Their  receiver  used  two  identical  units  and  a  common  local  oscillator  to  provide  time- 
aligned  phase  residuals.  This  is  essentially  the  same  approach  adopted  by  Ashtech  and 
various  independent  receiver  attitude  experiments. 

Adroit  -  Adroit  Systems  Inc.  has  developed  several  versions  of  an  attitude  determination 
system  for  artillery  alignment  [Jurgens  and  Rodgers,  1991]  and  for  land  and  marine 
vehicle  heading  determination  [Diefes,  1993].  The  original  design  was  a  set  of  three 
colinear  antennas  attached  to  a  single  receiver.  TTieir  newest  product  is  composed  of  four 
independent  receivers  which  collect  data  from  antennas  mounted  on  a  small  circular 
plate. 

Codeless  Attitude  Determination  System  -  Snyder  and  MacDoran  [1993]  at  the 
University  of  Colorado  have  designed  a  single  baseline  yaw  sensor  for  space  applications 
based  on  codeless  technology.  The  same  approach  could  be  extended  to  a  full  three  axis 
system.  The  technique  uses  separate  downconversion  stages  for  each  antenna  driven  by  a 
common  oscillator.  In  order  to  avoid  multipath  reflections  from  the  vehicle,  their  design 
incorporates  narrow  beam  antennas  (about  100  deg  field  of  view)  rather  than  the  typical 
hemispherical  patch  pattern. 

Trimble  -  The  Trimble  Vector  is  currently  the  only  receiver  that  employs  a  higher  level  of 
hardware  integration  than  other  approaches.  The  Trimble  Vector  is  a  six  channel,  C/A 
code,  single  frequency  receiver  which  multiplexes  among  four  antennas.  The  use  of  a 
single  hardware  processor  and  antenna  multiplexing  eliminates  any  changes  introduced 
by  independent  hardware.  Using  shared  hardware  for  all  RF  processing  adds  the 
complication  of  multiplexing  antennas.  This  introduces  difficulties  in  acquisition  and 
integer  tracking  and  reduces  the  effective  signal  to  noise  ratio  of  the  received  signal,  but 
these  issues  have  been  adequately  addressed  in  the  receiver  design. 

The  architecture  for  the  Vector  receiver  is  described  in  the  Trimble  Operations  Manual 
[Trimble,  1993]  and  Cohen  [1992].  The  top-level  design  for  the  receiver  is  shown  in 
Figure  5.1  [from  Cohen,  1992,  p.35].  The  RF  front  end  amplifies  the  signal  from  tlie 
antennas,  down-converts  the  signal  to  an  inteimediate  frequency  of  4  MHz  where  it  is 
digitally  sampled  and  passed  to  the  IQ  correlator.  Note  that  the  signals  from  all  four 
antennas  are  pass  through  the  same  RF  section  in  the  time-shared  multiplexing  scheme. 
From  here  the  digital  output  for  the  master  antenna  and  the  three  slave  antennas  are 
treated  differently.  The  master  antenna  output  drives  the  tracking  loop  filter  which 
regulates  the  numerically  controlled  oscillator.  All  signals  are  referenced  to  this  single 
local  oscillator,  and  the  oscillator  runs  open  loop  during  the  times  when  the  slave  antenna 
signals  are  being  tracked.  The  outputs  for  the  slave  antennas  are  passed  through  a 
separate  transfer  function  to  produce  the  final  differential  phase  output. 
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This  receiver  architecture  is  an  effective  signal  processing  design  for  measuring  the 
differential  phase  residual.  Nevertheless,  errors  are  still  present  in  the  raw  measurements. 
The  principal  sources  of  error  for  any  multi-antenna  GPS  receiver  are  receiver  biases, 
multipath,  antenna  cross-talk,  antenna  phase  center  motion.  Receiver  noise  includes 
random  measurement  errors,  tracking  errors,  antenna  cross  talk,  and  bias  drifts  due  to 
temperature  and  voltage  changes.  Multipath  and  antenna  phase  center  motion  will  be 
addressed  further  in  Section  5.6.  The  receiver  hardware  design  has  been  achieved  such 
that  the  receiver  errors  are  typically  less  than  2  mm  while  multipath  errors  can  be  as  large 
as  5  mm  [Braasch,  1992]. 

Trimble  also  manufactures  an  attitude  receiver  called  Quadrex  based  on  the  same  design 
as  the  Vector.  The  receivers  are  the  same  except  that  the  Vector  contains  an  additional 
internal  board  that  generates  attitude  solutions  without  external  processing.  The  Quadrex 
requires  an  external  computer  to  form  attitude  solutions. 

In  summary,  each  approach  has  its  own  challenges.  The  use  of  independent  receiver 
hardware  is  accomplished  most  efficiently  with  a  single  local  oscillator,  but  there  may  be 
additional  hardware  biases  introduced  by  temperature  or  voltage  changes.  The  common 
hardware  approach  adopted  by  Trimble  introduces  additional  complexities  in  the  signal 
processing  and  reduces  signal  to  noise  ratio.  Both  these  approaches  have  led  to 
successful  based  attitude  determination  systems  with  neither  proving  clearly 
superior  to  the  other. 


5.2  Least  Squares  Solution 

The  least  squares  algorithm  generates  a  discrete,  point  solution  (not  filtered)  for  the 
optimal  attitude  at  each  measurement  epoch.  The  Trimble  Vector  receiver  can  determine 
the  vehicle  attitude  using  a  nonlinear  least  squares  solution  or  the  SVD  solution  for  the 
Wahba  problem,  as  discussed  in  Section  4.6.  This  section  treats  only  the  nonlinear  least 
squares  solution  since  this  was  the  primary  method  used  in  subsequent  ground  tests. 
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The  attitude  solution  relies  on  additional  information  that  comes  from  an  off-line  baseline 
survey  and  resolution  of  the  integer  ambiguities.  Once  these  quantities  are  known,  the 
solution  is  found  from  the  linear  least  squares  equations, 

H50=5r  (5.1) 

where 


represents  small  angle  rotations  about  the  estimated  body  1,2,3  axes.  is  the  baseline 
cross  product  matrix  defined  in  Section  4.  The  attitude  matrix  estimate  can  be  updated 
using  50,  and  the  process  is  repeated  until  the  changes  get  sufficiently  small. 

Recall,  however,  that  the  GPS  attitude  receiver  measures  the  L-band  phase  difference 
between  the  master  and  each  of  the  baseline  antennas,  rather  than  the  range  difference.  If 
the  range  difference,  phase  difference,  and  baseline  vectors  are  all  described  in  cycles  (1 
cycle  -  19  cm  for  the  LI  frequency),  the  measured  phase  difference  is  just, 

Aq>ij=Arij-kij+Pi+Vij  (5.2) 

where  ky,  is  the  integer  ambiguity,  or  unknown  number  of  whole  cycles  in  the  phase 

difference.  Pi  is  a  line  bias  associated  with  the  receiver  hardware,  and  Vij  represents 
measurement  errors. 

Thus,  we  can  use  Equation  5.2  to  convert  a  GPS  phase  difference  measurement  to  a  range 
difference,  if  we  know  the  line  biases  for  each  slave  antenna  iPi,  for  i=l,..,/n)  and  the 
integer  ambiguity  for  each  satellite  and  slave  antenna  pair  (ity,  for  To  use 

the  range  differences  in  an  attitude  solution  we  also  require  the  baseline  vectors  in  the 
body  frame  (b/^,  for  i=l,../n  baselines),  and  the  line  of  sight  vectors  to  the  satellites  in 
the  navigation  frame  (ey^,  for  satellites).  The  body  frame  and  navigation  frame 

are  shown  in  Figure  5.2.  The  body  frame  has  the  x-axis  along  the  velocity  vector  of  the 
vehicle,  z-axis  pointing  up,  and  y-axis  perpendicular  to  both  in  a  right-handed  sense.  The 
navigation  frame  is  simply  the  earth-based  local  level  frame  expressed  in  East,  North,  Up 
components. 
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LOCAL  LEVEL  BODY 


Figure  5.2  Navigation  and  Body  Fixed  Reference  Frames  for  Ground  Tests 


Note  that  the  least  squares  cost  function  in  Equation  5.1  assumes  the  baseline  vectors  are 
known  and  that  differential  range  measurements  are  available.  The  next  section 
addresses  how  the  baseline  vectors  and  differential  range  values  are  derived.  If  the 
baseline  vectors  are  known  (and  assumed  fixed)  and  differential  range  measurements  are 
available  from  the  receiver  (i.e.,  the  integer  ambiguities  are  resolved),  then  only  two 
satellites  are  needed  to  form  attitude  solutions.  This  can  be  seen  by  examining  the  least 
squares  equation, 

HSe  =  Sr  (5.3) 

The  size  of  H  is  /nnx3,  Sr  is  mnxl,  and  the  size  of  Sd  is  3x1.  For  three  non-coplanar 
baselines  (m=3)  with  at  least  two  non-coincident  satellites  in  view  (n=2),  the  system  is 
overdetermined.  If  the  baselines  are  coplanar,  the  system  still  can  be  overdetermined  as 
long  as  the  satellites  do  not  like  in  the  plane  of  the  baselines. 


5.3  Baseline  Survey 

The  current  implementation  of  the  solution  algorithm  for  the  Trimble  receiver  requires 
three  distinct  steps:  1)  static  baseline  survey,  2)  integer  ambiguity  resolution,  and  3)  real 
time  attitude  solutions.  This  section  describes  a  static  baseline  survey  and  line  bias 
solution  technique. 

The  survey  is  required  to  define  the  baseline  vectors  in  body  coordinates  (b,®)  and  to 
determine  the  line  bias  for  each  slave  antenna  (^,).  The  geometrical  baseline  vectors 
could  be  measured  with  a  ruler;  however,  they  are  not  necessarily  coincident  with  the 
electrical  baselines.  In  addition  because  of  the  need  to  determine  the  line  bias,  the  GPS 
static  survey  is  a  convenient  means  for  establishing  both  the  vectors  and  line  biases.  To 
perform  a  survey,  the  vehicle  must  be  static  while  data  is  collected  from  several  satellites 
over  about  6-10  hours.  Long  satellite  arcs  are  required  to  separate  the  line  bias  from  the 
baseline  components. 
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If  the  satellite  line  of  sight  vectors  used  in  the  survey  are  expressed  in  East,  North,  Up 
(ENU)  components,  the  survey  will  solve  for  coordinates  of  each  slave  antenna  relative  to 
the  master  in  a  local  level  or  navigation  frame.  To  transform  them  to  the  body  frame,  the 
vehicle  must  be  aligned  in  a  known  (or  defined)  direction  relative  to  the  local  level.  The 
transformation  matrix  relating  the  two  frames  may  be  defined  by  calibration  or  by 
assuming  that  the  vehicle  is  aligned  in  a  given  direction. 

The  Vector  survey  software*  takes  the  latter  approach.  It  assumes  that  the  vehicle  is 
level,  i.e.  roll  and  pitch  are  zero,  with  baseline  2  serving  as  the  heading  reference.  Thus, 
once  the  baseline  vectors  are  determined  in  the  local  level  frame,  the  transformation  to 
body  coordinates  corresponds  to  a  single  rotation  (yaw)  about  the  vertical  to  align 
baseline  2  with  true  North. 

Each  master-slave  baseline  may  be  solved  independently.  The  number  of  unknowns 
which  must  be  determined  is  (N+3)  -  3  baseline  components  of  the  vector,  b„  and  one 
(ambiguity+line  bias)  for  each  of  the  satellites  tracked.  The  ambiguities  are  not  really  of 
interest,  but  the  common  fractional  part  determines  the  hardware  line  bias.  The  normal 
equations  are  formed  by  stacking  differential  phase  measurements  to  each  of  the  satellites 
over  time.  Of  course,  the  sequence  of  the  equations  does  not  matter  so  they  can  be 
grouped  conveniently  by  satellite  as  follows. 


Ax  =  5r 


where  for  n  satellites. 
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where  Ni  is  the  number  of  measurements  collected  for  satellite  /. 

Each  row  of  the  A  matrix  contains  a  line  of  sight  vector  to  one  of  the  satellites  and  a  "1" 
in  the  column  corresponding  the  integer+bias  for  the  particular  satellite.  The  hardware  or 
line  bias,  P  ,  is  defined  to  be  between  0  and  1,  and  is  taken  to  be  the  mean  (on  a  unit 
circle)  of  the  fractional  parts  of  (ki+p)  through  {kfi+P). 


*  The  jxototype  version  which  runs  on  an  external  computer  is  called  "SRVYTEST".  It  stores  a  file  called 
survey xam  which  contains  the  baslines  and  line  bias  values. 
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In  the  current  implementation,  several  hours  of  phase  measurements  are  typically 
processed  to  allow  enough  satellite  motion  to  reliably  estimate  the  fractional  hardware 
delay.  We  hope  that  future  implementations  of  this  algorithm  will  estimate  the  baseline 
vectors  and  possibly  the  line  biases  in  real  time,  eliminating  the  need  for  an  off-line  self 
survey.  This  would  be  valuable  for  spacecraft  applications  since  antenna  configurations 
computed  from  a  pre-launch  survey  could  change  slightly  due  to  stress  and  vibrations 
encountered  during  launch  and  deployment  Furthermore,  line  biases  may  drift  in  the 
space  environment  in  ways  not  accounted  for  in  a  pre-launch  survey.  Another  advantage 
of  estimating  the  baseline  vectors  in  real  time  is  that  baselines  might  deliberately  change 
if  antennas  are  placed  on  articulating  appendages  such  as  solar  panels.  In  this  case,  it 
may  be  possible  to  estimate  parameters  related  to  vehicle  configuration  as  well  as  vehicle 
attitude.  This  topic  is  highlighted  as  an  area  for  future  research. 

5.4  Integer  Ambiguity  Resolution 

Computing  the  integer  ambiguities  associated  with  each  antenna  and  satellite  is  a 
challenging  problem.  The  software  provided  with  the  Trimble  Vector  implements  a 
motion  based  approach,  as  described  in  Cohen  [1992]  and  Cohen  and  Parkinson  [1992a]. 
The  algorithm  for  integer  ambiguity  resolution  based  on  large  angle  motion  requires  the 
changes  in  the  baseline  vectors  be  accumulated  to  solve  for  the  baseline  vectors.  A  single 
baseline  change  for  an  arbitrary  number  of  baselines  is  shown  in  Figure  5.3.  The  large 
angle  motion  approach  is  based  on  the  observation  of  changes  in  the  differential  range 
between  antennas  expressed  as: 

(Arij  -  Ar,yo)  =  e /^x,-  (5.5) 

At  each  epoch  Ax,  must  be  solved  for  each  of  the  m  baselines  by  forming  a  solution 
using  a  minimum  of  three  satellites  in  view  and  assuming  the  line  of  sight  vectors  are 
constant  for  the  short  time  period  considered.  This  can  be  expressed  in  matrix  form  as: 

Axj.  (5.6) 


Fig.  5.3.  Baseline  Changes  for  Ambiguity  Resolution 
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Cohen  [1992]  points  out  that  the  constraint  on  the  entire  baseline  structure  can  be  used  to 
resolve  the  integer  ambiguities.  Since  the  baselines  are  fixed  (both  in  length  and  relative 
orientation), 

Xj(ti)  =  XiUo)^  Xjito) 

for  each  pair  of  i  and  j. 


If  we  express  baseline  i  at  rj  as  =  x,o  +  Ax,  (ri) .  then 

X,  iOf  Xj  (0)  =  [x,  (0)  +  At,  (/i  )]^[x^  (0)  + Ax  j  (fi )] 


We  need  to  find  x,(ro)  and  xjitQ)  in  order  to  solve  for  the  integer  ambiguities.  We  have 
computed  Ac,  from  Equation  5.5,  so  we  can  combine  N  measurements  from  20-60 
epochs  in  matrix  form  to  solve 

Ax/»)^  1  r  Ax,„(l)^  1  r  Ac/1)7-At„«> 

^m+  :  Jf/=- 

Ac/^)7-J  [ax„<^)^J  [ax/^>^Ac„(^) 

where  superscripts  in  parenthesis  denote  each  epoch.  These  equations  can  be  combined 
into  a  single  linear  least  squares  equation  of  the  form  Ax^B  to  solve  for  the  baseline 
vectors,  x.  The  specific  formulation  of  this  single  equation  depends  on  the  number  of 
baselines  used. 

Once  the  initial  baseline  vectors  are  known  each  of  the  kjj  terms  can  be  found  using 
and  the  precomputed  Pi  from  Equation  5.2. 


(5.9) 


A<Pij^Ar.j-kij+Pi 

where  the  noise  term  Vij  has  been  assumed  to  be  zero. 

The  ambiguities  are  then  rounded  to  integers  and  the  standard  solution  algorithm  can 
begin. 

The  Trimble  Vector  receiver  uses  the  large  angle  integer  resolution  approach  based  on  the 
information  from  the  static  self  survey  to  initiate  the  attitude  solution.  During  normal 
operation,  the  integers  are  regularly  recomputed  to  enable  new  satellites  to  be  brought 
into  the  solution  and  to  verify  the  current  integer  values.  If  the  residuals  are  found  to  be 
too  large,  the  algorithm  will  "bump"  them  to  reduce  the  residuals  to  an  acceptable  level. 
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5.5  (iround  Testing  Configuration 

Ground  lusung  ot’  the  Tnmrie  X'cetor  receiver  '.vas  u.M-ieuctee  enuer  retn  :aue  ana 
dynamic  i.^pinnlngi  tundiUun.s.  The  tests  toeused  on  measurement  uuaiitv,  static  sen 
>urve>  .  static  attitude  solutu  :■  accuracs,  ana  dsnamic  reriormancc  .  n  ..  spinning 
piattorm.  Before  addressing  the  specific  experiments,  tne  test  Peas  ssiii  re  aesenneu. 
Each  test  bed  consisted  of  a  site,  an  antenna  platform,  ^nd  known  reference  points  to 
verify  accuracy. 

Test  locations 

Two  locations  were  u.sed  lor  gruunu  testing  --the  roof  ot  the  Engineenng  Center,  and 
Table  Mountain.  The  roof  of  the  Engineenng  Center  at  the  Cniversity  of  Colorado 
eonsists  of  a  flat  area  approximately  60  x  100  m  (.see  Figure  5.4).  The  view  of  the  sky  is 
partially  obstructed  to  the  southwest  by  the  tower  of  the  Engineenng  Center  and  adjacent 
buildings.  Many  antennas  and  metal  structures  are  co-loeated  on  the  roof,  making  it  a 
very  noisy  environment  for  multipath.  (However,  it  was  quite  convenient. ) 


Fig.  5.4.  Engineering  Center  Roof  and  x  Structure 

The  other  test  location  used  was  Table  .Mountain,  shown  in  Figure  5.5.  Table  Mountain 
is  a  large  flat  mesa  north  of  Boulder.  Colorado,  controlled  by  the  U.S.  Department  of 
Commerce.  The  area  used  was  largely  free  of  multipath  reflectors  and  contained 
surveyed  reference  points. 
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Fiu.  5.^  Table  Mouniain  Sue  aaia  Tree  Suiieture 


Antenna  Platforms 

Four  antenna  platforms  were  used  for  receiver  lesiing.  Each  platform  was  tailored  for  a 
particular  purpose  in  the  testing  program.  The  platforms  will  be  referred  to  as  the  x 
structure,  the  tree  structure,  the  one-armed  spinner,  and  the  wooden  palleLs. 

The  X  structure  is  constructed  from  two  aluminum  cross  members,  each  1-m  long  as 
.shown  in  Figure  5.4.  The  master  antenna  is  raised  on  spacers  that  are  approximately 
10  cm  long  to  prevent  all  four  antennas  from  being  coplanar.  .A  small  diode  laser  can  be 
mounted  at  the  center  of  the  two  cross  members  for  absolute  accuracy  tests.  The  x 
structure  was  used  in  early  tests,  but  it  suffered  from  several  shoncomings  --the  ba.selines 
were  too  short;  the  antennas  were  too  clo.se  to  ci  "“lanar:  the  arms  were  too  flexible;  and 
the  structure  itself  was  a  multipath  reilector. 

The  tree  structure  shown  in  Figure  5.5  was  designed  to  address  the  problems  with  the  x 
structure.  It  has  a  center  mounting  pole  for  the  master  antenna  and  three  one-meter  arms 
for  the  slave  antennas.  The  center  pole  and  baseline  2  define  the  body  coordinate  system. 
The  arms  can  be  moved  to  change  the  relative  height  and  orientation  of  the  slave 
antennas.  The  structure  is  fabricated  using  Schedule  40  PVC  pipe.  This  material  was 
cho.sen  because  it  was  thought  to  have  several  advantages  (Comp,  et  al..  IdO.’^j  -  greater 
rigidity,  preventing  significant  flexure  during  the  tests;  lighter  weight,  permitting  e.^sy 
relocation  and  storage;  and  transparency  to  L-band  GPS  signals.  The  small  diode  laser 
was  mounted  below  the  baseline  2  antenna  for  absolute  accuracy  tests.  Note  that  the 
three  antenna  baselines  shown  in  the  photograph  are  nearly  orthogonal  with  the  master 
antenna  at  the  vertex  of  the  orthogonal  set.  This  nearly  orthogonal  configuration  provides 
the  best  sensitivity  to  changes  in  attitude  for  all  three  axes.  Field  experience  showed  that 
the  tree  structure  was  not  as  rigid  as  originally  hoped,  thus  suggesting  the  area  of  attitude 
determination  for  non-rigid  structures  as  an  area  for  future  research. 
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The  v)no-armcd  spinner  Niruciure  was  aesiiznea  ui  nreuni  vP.  ^n  exiswiie  -pinnir.e  .:v 

The  siruciurc  consists  •  !  a  sineie  mastcr-siave  antenna  r.e.r  rr.eiuniee  en  a  P\X'  .ante; 
pipe  and  a  sineie  i  meter  ione  PV’C  pipe  arm.  The  arm  uas  t'litea  to  a  e>iinu^r  .maante. 
-  n  the  spin  table,  Fne  lanle  anven  ny  a  DC  inoaor  -a  nun  :'.a.s  a  .iriaPiL-  ppAver  .'■appi' 
:hat  can  he  used  to  reeuiate  spin  rate.  The  receiver,  pov^er  suppiv.  ana  Jau  ilee:;.  ; 
computer  are  all  stored  inside  i.Pe  cylinder  aunne  lesime.  "he  eoaii leuration  ;>  snoun  ir 
Fieum  5.P. 


Fig.  5.6  One-Armed  Spinner  Structure 
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The  wooden  palleus  were  usc'd  lor  tests  wnere  antenna  loeatton  ana  onentation  naa  to  he 
maintained  preeisely  t’rom  one  dav  tu  the  next.  Twi>  wiu'den  pallets  sveie  plaeea 
approximately  3  m  apart  on  the  Hat  root'  oi'  a  low  huildine  at  Table  .Mountain.  This 
eliminated  any  structural  tlexure  that  mav  have  been  pre.vmi  in  the  otter  '.e.^t  .strueture.s. 
The  master  antenna  was  m.'unted  to  one  oi'  the  pallets  with  C-eiamp.s.  and  two  slave 
antennas  were  mounted  to  the  other  pallet.  Fiszure  y.”  shows  the  pallet  eonfi duration  at 
Table  Mountain. 


Fig.  5.7  Wooden  Pallet  Structure 
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Reference  Points 

Known  control  points  at  the  Table  Mountain  test  site  provided  a  ground  truth  reference 
for  attitude  tests.  Table  Mountain  has  been  used  by  UNAVCO  for  training  personnel  in 
GPS  surveying  techniques  [Roberts,  1990;  UNAVCO,  1993].  The  vectors  connecting  the 
control  points  establish  known  directions  in  an  Earth  fixed  or  local  level  coordinate 
system.  Absolute  attitude  determination  accuracy  using  GPS  may  be  evaluated  by 
centering  the  structure  on  one  control  point  and  aligning  it  with  another.  For  the  relative 
and  absolute  accuracy  tests  described  below,  the  antenna  structure  was  placed  at  the  base 
point  shown  in  Figure  5.8.  The  coordinates  of  the  other  points  relative  to  the  base  point 
are  shown  in  Table  5.1. 


Fig.  5.8  Table  Mountain  Control  Points 

Table  5.1-  Relative  Positions  of  Control  Points  from  Mount  (from  [UNAVCO,  1993]). 


SITE 

EAST 

(m) 

NORTH 

(m) 

UP 

(m) 

AZ 

(dCR) 

EL 

(deg) 

RHO 

(m) 

A 

35.667 

36.451 

-0.544 

44.376 

-0.612 

51.001 

B 

^1.922 

55.465 

1.029 

322.917 

0.842 

69.533 

C 

-3.551 

34.879 

0.273 

354.186 

0.446 

35.061 
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Basic  Operations 

To  initialize  attitude  operation,  the  equipment  is  set  up  as  shown  in  Figure  5.9.  The 
equipment  consists  of  the  Trimble  receiver,  power  supply,  computer,  and  an  appropriate 
antenna  mounting  structure.  Receiver  initialization  commands  are  issued  from  the 
computer.  The  structure  is  then  rotated  by  about  90  degrees  to  resolve  the  integer 
ambiguities,  and  the  attitude  solutions  are  then  initialized.  Data  files  are  stored 
containing  the  raw  differential  phase  information,  direction  cosines  for  each  of  the 
tracked  satellites,  and  attitude,  position,  and  velocity  solutions.  The  time  interval  for 
storage  is  internally  set  to  0.5  seconds  by  the  receiver,  but  this  value  may  be  adjusted  by 
the  user. 


Fig.  5.9.  GPS  Receiver  and  Test  Structure  Configuration 
5.6  Phase  Measurement  Quality 

The  characteristics  of  the  phase  residuals  were  examined  to  determine  the  quality  of  the 
raw  measurements.  Multipath  and  antenna  phase  center  migration,  potentially  the  two 
largest  error  sources,  were  considered  separately  in  assessing  overall  measurement 
quality.  Multipath  effects  were  studied  through  both  experiments  and  modeling. 

The  quality  of  the  raw  differential  phase  observations  can  be  assessed  by  removing  the 
effect  of  satellite  motion  leaving  only  the  "phase  residual".  If  the  location  and  orientation 
of  the  baseline  is  known,  the  effect  vif  the  satellite  motion  can  be  computed  based  on  this 
information  and  the  transmitted  satellite  ephemerides.  This  predicted  phase  difference  is 
subtracted  from  the  observed  phase  difference  to  form  the  residual.  If  the  orientation  is 
not  known,  a  polynomial  of  second  or  third  order  (depending  on  the  length  of  the  data 
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set)  may  be  removed  to  detrend  the  data.  The  latter  method  was  used  in  most  of  the 
results  presented  here. 

Multipath  Testing 

Figure  5.10  shows  the  satellites  in  view  for  a  data  set  collected  on  May  9,  1993.  The  test 
was  conducted  at  Table  Mountain  using  the  tree  structu.<e.  Figure  5.11  illustrates 
measurement  residuals  for  GPS  satellite  2  from  each  of  the  three  baselines.  As  shown  in 
the  figure,  the  random,  or  noise-like  error  in  the  observed  residuals  was  on  the  order  of  2 
mm,  1-0.  Over  the  1.4  meter  baseline  used  in  the  tree  structure,  this  would  produce  a 
pointing  error  of  less  than  0.1  degrees.  These  extremely  good  results  can  be  attributed  to 
the  high  elevation  (shown  in  Figure  5.10)  for  this  satellite. 


SATELLITES  IN  VIEW  DURING  RUN  M09TR0N2 


Fig.  5.10.  Azimuth-Elevation  Plot  of  Satellites  in  View  for  Test  at  Table  Mountain,  9 
May  1993.  Radial  distance  represents  satellite  zenith  angle  (90  deg  -  elevation).  Shaded 
hemisphere  indicates  the  position  of  a  large  metal  dome. 

Systematic  errors  can  also  be  observed  in  some  of  the  phase  residual  data.  Figure  5.12 
illustrates  the  results  for  satellite  13  during  the  same  May  9,  1993  test.  From  Figure  5.10 
it  can  be  seen  that  satellite  13  was  at  an  elevation  angle  of  only  about  15-20  degrees.  In 
addition,  it  was  in  the  direction  of  a  large  metal  dome  at  the  test  site.  Systematic  errors 
such  as  these  are  probably  due  to  multipath  reflection  of  satellite  signals  from  nearby 
objects  or  variations  in  the  phase  center  of  the  microstrip  antennas. 
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09  MAY  93  RUN#2  ROOF  ->  (NORMAL.TBL  )  SV2 


Fig.  5.11.  Phase  Difference  Residuals  for  PRN  2,  9  May,  1993.  Solid  line  shows 
baseline  1  residuals  offset  by  0.15  cycles;  dotted  line  shows  baseline  2  residuals;  dashed 
line  shows  baseline  3  residuals  offset  by  -0.15  cycles. 

09MAY93RUN#2  ROOF“>  (NORMALTBL )  SV13 
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Fig.  5.12.  Phase  Difference  Residuals  for  PRN  13,  9  May,  1993.  Solid  line  shows 
baseline  1  residuals  offset  by  0.15  cycles;  dotted  line  shows  baseline  2  residuals;  dashed 
line  shows  baseline  3  residuab  offset  by  -0.15  cycles. 
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Systematic  errors  may  be  repeatable  for  certain  aerospace  applications.  For  example,  a 
satellite  in  orbit  may  closely  repeat  its  orientation  with  respect  to  the  GPS  constellation. 
Several  tests  were  conducted  to  see  if  the  repeatable  multipath  similar  to  those  reported 
by  Cohen  and  Parkinson  [1991a]  and  Braasch  [1992]  could  be  observed  in  the  raw  phase 
data. 

Tests  for  multipath  repeatability  were  conducted  on  the  roof  of  the  Engineering  Center 
using  the  x  structure.  Phase  residuals  were  collected  on  several  dates  with  5  I^b  93  taken 
as  a  reference.  Data  from  the  following  days  were  shifted  by  4  minutes  per  day  to 
compensate  for  the  advance  of  the  satellite  ground  tracks.  Line  of  sight  vectors  to  the 
satellites  in  view  were  compared  to  confirm  that  the  same  satellite  alignment  was 
achieved  for  the  various  data  sets.  Figure  5.13  shows  phase  residuals  for  SV  15  collected 
two  weeks  apart.  Satellite  motion  has  been  removed  from  the  residuals  by  a  quadratic  fit 
and  the  data  for  the  three  baselines  was  offset  to  display  on  the  same  graph.  The  satellite 
was  fairly  low  on  the  horizon  during  these  tests  (~  25  degrees),  and  it  was  setting  behind 
the  tower  of  the  Engineering  Center.  Comparison  of  the  results  from  5  Feb  and  19  Feb 
93  shows  a  very  strong  repeatability  of  the  multipath  effect  for  all  three  antenna 
baselines. 


05  FEB  93  RUN#1  MOUNT->  (NORMAL.ROOF)  SV15 
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19FEB93  RUN#1  MOUNT->  (NORMAL.ROOF)  SV15 
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TIME  (min) 

Fig.  5.13.  Phase  Difference  Residuals  for  PRN  15.  (a)  is  for  5  Feb  1993,  (b)  is  for  19 
Feb  1993.  Solid  line  shows  baseline  1  residuals  offset  by  0.15  cycles;  dotted  line  shows 
baseline  2  residuals;  dashed  line  shows  baseline  3  residuals  offset  by  -0.15  cycles. 
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Repeatability  of  multipath  effects  for  these  tests  was  also  analyzed  in  the  frequency 
domain  by  examining  the  power  spectrum  of  each  of  the  runs.  The  raw  phase  data  from 
two  different  days  were  aligned  to  have  the  same  satellite  geometry,  as  before.  A  Fast 
Fourier  Transform  (FFT)  was  used  to  generate  power  spectra  for  the  two  cases.  Figure 
5. 14  shows  the  power  spectral  densities  for  30-minute  runs  on  5  and  19  Feb.  Frequencies 
are  normalized  to  the  sampling  rate  of  2  Hz.  Energies  at  specific  frequencies  are  not 
evident  from  these  plots  and  there  is  no  obvious  correlation  between  the  two  dates.  This 
lack  of  clearly  defined  frequency  content  is  probably  due  to  Doppler  broadening  of  the 
characteristic  multipath  frequencies  due  to  satellite  motion.  Analysis  of  shorter  time 
segments  was  then  used  to  minimize  the  effect  of  satellite  motion. 


Pxx  -  X  Power  Spectral  Density 


Fig.  5.14.  Power  Spectrum  for  30-minute  Run  on  19  Feb  1993.  (a)  is  for  5  Feb  1993,  (b) 
is  for  19  Feb  1993.  Frequency  is  normalized  to  the  sample  rate  of  2  Hz. 
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Figure  5.15  shows  power  spectra  for  a  shorter  (5  min.)  portion  of  the  same  data  set.  Note 
that  there  are  more  "peaks"  in  the  graphs,  indicating  energies  at  discrete  frequencies. 
Also,  the  arrows  show  peaks  that  seem  to  be  repeated  for  the  two  runs,  possibly 
indicating  a  frequency  contribution  from  a  specific  multipath  reflector.  This  indicates 
some  repeatable  frequency  components  of  the  multipath  effect  which  could  be  due  to 
specific  combination  of  multipath  reflectors,  however,  the  results  are  certaintly  not 
conclusive. 


Pxx  -  X  Power  Spectral  Density 


Fig.  5.15.  Power  Spectrum  for  5-Minute  Run  on  5  Feb  1993.  (a)  is  for  5  Feb  1993,  (b)  is 
for  19  Feb  1993.  Frequency  is  normalized  to  the  sample  rate  of  2  Hz. 
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Multipath  Modeling 

Analytical  modeling  of  multipath  was  used  to  try  to  identify  the  effects  of  single 
reflectors  of  the  GPS  signal  based  on  geometric  optics.  Geometric  optics  consists  of  ray 
tracing  for  each  reflector  in  the  multipath  environment. 

Expressions  for  the  composite  carrier  phase  error  due  to  a  single  multipath  reflector  and 
a  fixed  satellite  geometry  have  been  derived  in  Braasch  [1992].  These  expressions 
essentially  account  for  a  fixed  frequency  disturbance  in  the  carrier  phase  that  depends  on 
the  distance  from  the  reflector  to  the  receiver  and  the  relative  magnitudes  of  the  direct 
and  reflected  signals.  Figure  5.16  shows  the  power  spectrum  for  simulated  phase 
residuals  from  a  single  multipath  reflector  using  the  expressions  found  in  Braasch  [1992]. 
The  elbow  in  the  spectrum  at  a  (normalized)  frequency  of  about  0.08  is  the  energy  due  to 
the  multipath  reflector. 


Pxx  -  X  Power  Spectral  Density 


Fig.  5.16.  Power  Spectrum  for  Simulated  Phase  Residuals  Using  Static  Model 

The  difficulty  with  this  model  is  that  it  does  not  account  for  changes  in  the  path 
difference  due  to  satellite  motion.  An  analytical  model  for  a  single  multipath  reflector 
was  developed  that  does  account  for  changing  geometry  due  to  satellite  motion.  The 
basic  geometry  is  shown  in  Figure  5.17.  The  equations  for  path  length  difference 
between  the  direct  and  reflected  signals  for  a  reflector  placed  a  distance  d  from  the 
antenna  can  be  shown  to  be: 
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6=  V-  sin* 


^sin(^-u) 


'  =  -^raJ^  +  re^-2rat,  '■gCOS 


6  =  Vr^  +  -2rdco8  0 


AX  =  b  +  d-r  (5.10) 

where,  AX  is  the  path  length  difference,  r  is  the  path  length  for  the  direct  signal,  is  the 
radius  of  the  earth,  is  the  satellite  orbit  radius,  6  is  the  satellite  look  angle,  and  v  is  the 
satellite  mean  motion.  All  lengths  are  expressed  in  LI  carrier  wavelengths. 


Figure  5.17  -  Dynamic  Multipath  Model  Geometry 

The  time  difference  between  the  arrival  of  the  direct  and  reflected  signal  varies  as  the 
spacecraft  travels  overhead.  The  delay  attains  a  minimum  (zero)  when  the  reflector  is 
between  the  satellite  and  the  receiver  and  all  three  are  collinear,  and  reaches  a  maximum 
when  all  three  are  collinear  with  the  satellite  and  receiver  on  opposite  sides  of  the  receive 
antenna. 

Simulated  multipath  errors  based  on  this  model  were  compared  with  experimentally 
obtained  phase  residuals.  These  attempts  to  correlate  ray  tracing  techniques  with 
experimental  results  proved  frustrating.  Experimental  results  represent  a  composite  of  all 
multipath  reflectors  in  the  local  environment,  making  it  nearly  impossible  to  isolate  the 
effects  of  a  single  reflector  (or  even  a  few  reflectors).  Thus,  although  geometric  optics  is 
a  useful  tool  for  illustrating  the  basic  mechanism  of  multipath,  we  were  not  able  to  use  it 
to  evaluate  experimental  results. 

Another  way  of  gaining  insight  into  this  effect  that  has  been  proposed  by  Cohen  [1992, 
and  Cohen  and  Parkinson,  1991a]  is  a  spherical  harmonic  expansion  of  the  multipath 
environment  derived  from  empirical  results.  Although  this  appears  to  work  well  under 
some  test  conditions,  it  requires  that  multipath  be  exactly  repeatable  for  a  given  satellite 
geometry.  This  may  not  be  realistic  if  multipath  reflectors  change  their  reflective 
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characteristics  under  different  temperature  or  moisture  conditions  or  they  change  their 
configuration  with  respect  to  the  antennas,  as  seen  in  the  next  experiment. 

Antenna  Phase  Center  Migration 

Another  source  of  systematic  error  that  must  be  considered  is  the  migration  of  the 
electrical  center  of  the  antenna  depending  on  the  direction  of  the  incoming  signal.  This 
effect  could  introduce  a  slowly  varying  error  in  the  results  that  depends  on  the  properties 
of  the  particular  antennas  used.  The  magnitudes  reported  for  these  errors  range  from 
about  0.8  cm  [Lightsey,  et  al.,  1993]  to  as  much  as  2  cm  [Davis,  et  al.,  1993].  A  test 
designed  to  isolate  antenna  phase  center  motion  for  our  antennas  was  conducted  at  Table 
Mountain  using  the  wooden  pallet  antenna  platform. 

The  basic  approach  for  isolating  antenna  phase  center  motion  was  to  design  a  field  test 
where  all  conditions  would  repeat  from  day  to  day  with  the  exception  of  antenna 
orientation.  By  comparing  these  results,  all  effects  other  than  antenna  phase  center 
motion  would  cancel  from  the  collected  phase  residuals.  The  tests  were  conducted  with 
one  baseline  (Red-Blue)  fixed  as  a  control  for  all  four  runs  of  the  experiment  on  13-16 
July  1993.  The  antenna  for  the  other  baseline  (Red-Green)  was  rotated  on  successive 
days  as  shown  in  Figure  5.18. 


JULY  13  JULY  14 


Fig.  5.18.  Antenna  Configurations  for  Phase  Center  Migration  Tests.  Small  circle 
indicates  position  of  the  cable  connector.  The  red  "R"  antenna  is  the  master;  the  blue 
"B"  antenna  is  the  control  slave;  the  green  "G"  antenna  is  the  slave  antenna  under  test. 

The  recorded  phase  residuals  were  shifted  by  four  minutes  per  day  to  account  for  satellite 
motion.  The  data  was  then  low-pass  filtered  to  reduce  the  high  frequency  receiver  noise. 
Satellite  line  of  sight  vectors  stored  every  30  seconds  were  interpolated  to  the  one  second 
measurement  epochs.  The  phase  residuals  from  July  16, 1993  were  taken  as  the  reference 
data  set,  and  these  were  differenced  with  the  phase  residuals  from  the  other  days. 

Unfortunately  the  results  were  not  very  conclusive  because  of  the  difficulty  in  achieving 
sub-centimeter  repeatability  of  multipath  effects  in  a  less  than  ideal  outdoor  environment. 
For  example,  multipath  reflectivity  of  the  ground  and  other  objects  can  change  due  to 
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differing  moisture  levels.  Thus,  it  is  possible  that  the  results  from  July  13  and  14  were 
substantially  different  from  July  IS  and  16  because  of  rain  showers  on  the  earlier  two 
days  and  dry  weather  during  the  latter. 

Figure  5.19  shows  the  phase  residuals  for  July  16  subtracted  from  July  IS.  The  top  graph 
(Figure  S.  19(a))  shows  the  baseline  with  the  rotated  antenna  configuration  on  the  two 
days.  Figure  S.  19  (b)  shows  the  results  for  the  baselines  with  antennas  fixed  in  the  same 
position  on  both  days.  Note  that  the  top  graph  for  the  rotated  antenna  baseline  shows 
more  systematic  features  and  apparently  larger  deviations  from  zero.  The  standard 
deviation  for  the  rotated  antenna  baseline  is  3.7  mm  while  the  standard  deviation  for  the 
fixed  baseline  is  2.9  mm. 
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Fig.  5.19.  Antenna  Phase  Center  Migration  Test  Results  for  PRN  2  Comparing  15  Jul 
1993  and  16  Jul  1993.  (a)  is  for  baseline  with  rotated  antenna  configuration,  (b)  is  for 
baseline  with  fixed  antenna  configuration. 

For  all  satellites  in  view  during  this  test,  the  difference  in  standard  deviation  between  the 
rotated  and  the  fixed  antenna  baselines  ranged  from  0.1  mm  to  1.6  mm.  The  satellites  in 
view  during  the  test  are  shown  in  Figure  5.20.  Typically,  the  rotated  antenna  baseline 
standard  deviation  is  about  0.6  -  0.7  mm  larger  than  the  fixed  baseline.  The  smaller 
differences  would  indicate  a  small  phase  center  migration  for  the  particular  satellite 
motion  encountered,  while  the  larger  differences  indicate  larger  migrations  for  that 
combination  of  antenna  rotation  and  satellite  motion. 


PHASE  CENTER  MIGRATION-FIXED  BASELINE:  sv2 
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SATELLITES  IN  VIEW  DURING  RUN 


,100  L - - - - - i 

5.20.  Azimuth-Elevation  Plot  Satellites  in  View  During  Ar.tenra  Phase  Center  Motion 
Tests  at  Table  Mountain.  Radial  distance  represents  satellite  zenith  angle  (90  deg  - 
elevation). 

The  important  result  from  this  test  is  that  the  rotated  baseline  exhibited  larger  standard 
deviations,  and  the  graphs  typically  show  more  systematic  effects  than  the  fixed  baseline. 
This  strongly  suggests  that  antenna  phase  center  migration  is  being  detected  with  this  test 
Note,  however,  that  antenna  phase  center  migration  is  not  generally  as  large  as  multipath 
effects.  These  results  also  show  that  antenna  phase  center  migration  for  our  equipment  is 
much  smaller  than  the  20  mm  for  the  test  configuration  used  in  Davis,  et  al.  [1993]. 
Further  tests  could  be  used  to  calibrate  antenna  phase  center  migration  for  extremely 
accurate  attitude  solution  processing,  but  this  should  not  be  a  top  priority  until  multipath 
errors  have  been  reduced  through  some  means  to  the  same  order  of  magnitude  as  the 
antenna  phase  center  migration. 

5.7  Baseline  Survey  Test 

Static  tests  of  the  off-line  survey  were  conducted  to  verify  repeatability  of  baseline  vector 
and  line  bias  calculations.  The  goal  of  these  tests  was  to  identify  the  issues  that  will  need 
to  be  resolved  for  real  time  implementation  of  the  baseline  survey. 
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Equipment  and  Procedure 

Testing  was  conducted  using  the  wooden  pallet  structure  on  the  roof  of  the  engineering 
center.  The  pallets  ensured  that  structural  vibration  and  antenna  motion  could  be 
eliminated  during  the  baseline  surveys.  Two  separate  surveys  of  a  single  antenna 
baseline,  were  conducted  on  27  Sep  93.  They  will  be  referred  to  as  Survey  1  and  Survey 
2,  and  they  each  lasted  more  than  6  hours.  The  baseline  and  antenna  orientation  remained 
fixed  throughout  both  surveys.  The  data  was  post  processed  using  the  algorithm 
discussed  in  Section  5.3  to  provide  an  initial  estimate  of  the  baseline.  Data  outliers  were 
then  discarded  by  reprocessing  the  data  using  only  points  with  differential  range  residuals 
less  than  2  cm.  The  line  bias  for  the  baseline  was  then  determined  using  the  data  set  with 
outliers  removed. 

Results 

Table  5.2  shows  the  least  squares  baseline  solutions  for  the  two  surveys.  (Note  that 
measurement  outliers  were  removed.)  The  agreement  in  each  component  between  the 
runs  is  6  mm  or  better.  The  baseline  lengths  agree  to  better  than  2  mm. 


Table  5.2.  Baseline  Results  for  Survey  1  and  Survey  2 


Survey  1 

Survey  2 

East  Component 

4.0189  m 

4.0194  m 

North  Component 

1.2837  m 

1.2820  m 

Up  Component 

-1.3748  m 

-1.3800  m 

Length 

4.4373  m 

4.4389 

Table  5.3  shows  the  least  squares  solutions  for  the  integer  +  line  bias  values  for  each 
satellite  used  in  the  solution.  These  results  were  also  calculated  with  anomalous  data 
points  removed.  The  integers  differ  from  survey  to  survey,  as  expected,  since  their  value 
is  run  dependent.  However,  the  line  bias  calibration  should  remain  constant  between 
surveys.  Note  from  Table  5.3  that  six  satellites  were  used  in  Survey  1  while  only  four 
were  used  for  Survey  2. 


Table  5.3.  Integer  Plus  Line  Bias  Results  for  Survey  1  and  Survey  2.  Integer  plus 
_ bias  values  are  given  in  units  of  LI  carrier  wavelengths. _ 


Survey  1 

Survey  2 

16.895 

i099 

51^5 

0.877 

5095 

ros5 

5r9oo 

The  average  line  bias  for  survey  1  was  0.912  cycles,  and  the  average  for  survey  2  was 
0.890  cycles. 
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The  surveys  showed  sub-centimeter  repeatability  in  each  of  the  baseline  vector 
component  directions,  and  the  average  line  biases  agree  within  about  0.02  cycles  (less 
than  4  mm).  This  excellent  repeatability  provides  confidence  in  the  results;  however  the 
survey  algorithm  iS  not  completely  automated.  To  get  reliable  results  we  were  required  to 
manually  identify  outliers  and  eliminate  satellites  which  were  only  tracked  for  a  short  arc. 
This  does  not  present  a  problem  for  off-line,  post-processed  solutions,  but  treating  such 
issues  on-line  will  be  more  difficult 

Algorithms  for  real  time  implementation  of  the  baseline  survey  must  overcome  these 
challenges.  Anomalous  da.a  points  must  be  discarded  as  they  are  collected,  and  cycle 
slips  corrected  as  they  occur.  Since  the  duration  of  the  runs  is  long  (over  5  hoors)  many 
satellites  ^nter  and  leave  the  field  of  view  during  the  survey.  Satellites  which  are  only 
tracked  for  a  short  arc  cannot  be  used  to  determine  the  line  bias.  For  example,  consider 
the  last  row  in  Table  5.3  for  survey  1.  This  integer  plus  line  bias  value  corresponds  to  a 
satellite  that  was  rising  less  than  2  hours  before  the  survey  was  terminated,  resulting  in 
fewer  data  points  being  available  to  generate  this  value.  If  this  value  is  discarded  from 
the  average,  then  the  average  line  bias  for  sup/ey  1  is  0.892  cycles,  a  difference  of  only 
0.002  cycles  from  survey  2. 

On  the  other  hand,  discarding  too  many  data  points  or  eliminating  satellites  from  the 
algorithm  can  quickly  lead  to  underdetermined  solutions  for  the  baselines  and  integer 
plus  bias  values.  Automating  these  trade-offs  in  a  robust  real-time  algorithm  will  be  a 
major  challenge  in  future  research. 


5.8  Relative  Pitch  Accuracy  Test 

The  next  experiment  performed  was  designed  to  compare  the  change  in  pitch  reporuid  by 
the  receiver  to  an  independent  measure  of  the  angular  change  determined  by  the  motion 
of  the  laser. 

Equipment  and  Procedure 

The  tree  structure  is  fitted  with  the  laser  mounted  below  the  baseline  2  slave  antenna. 
The  laser  is  initially  pointed  at  the  wall  of  a  small  storage  shed  approximately  50  meters 
away.  The  laser  image  on  the  wall  (approximately  4  cm  in  diameter)  is  marked  and  the 
GPS  based  attitude  solutions  are  collected  for  about  15  minutes.  The  structure  is  then 
pitched  upward  until  the  laser  is  aligned  with  a  point  higher  on  the  wall.  The  higher  spot 
is  marked,  a.nd  data  is  again  collected  for  about  15  minutes.  Finally,  the  structure  is 
removed  from  the  control  point  and  replaced  with  a  Wild  T2  Theodolite  to  perform  an 
independent  measure  of  the  angular  distance  between  the  center  of  the  marks.  The 
estimated  precision  of  the  laser  image  marks  is  about  1  cm,  producing  an  angular  error  of 
0.0  i 2  degrees  on  a  50  meter  baseline.  The  precision  of  the  theodolite  measurement  is 
better  than  0.001  degrees.  Thus,  the  pitch  reference  angle  is  good  to  about  0.012  degrees. 
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Results 

Figure  5.21  shows  the  GPS  based  relative  pitch  angle  estimates.  The  pitch  values  are 
offset  so  that  the  mean  for  the  first  15  minutes  is  zero  degrees. 

The  measured  relative  change  in  pitch  during  the  test  was  2.07  deg.  Using  a  Wild  T2 
Theodolite  the  angular  separation  of  the  marics  was  measured  at  1.96  degrees,  indicating 
a  relative  accuracy  of  approximately  0.11  deg  for  this  test.  The  T2  Theodolite  optical 
measurement  is  accurate  to  better  than  0.001  deg. 

This  result  is  consistent  with  our  expectations  based  on  the  measurement  quality  and  the 
baseline  length. 
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5.21.  Real  Time  Relative  Pitch  Estimates,  10  May  1993 
5.9  Absolute  Heading  Accuracy  Test 

The  primary  goal  of  the  performance  tests  was  to  determine  the  accuracy  of  the  real  time 
attitude  solutions  as  compared  to  the  known  pointing  direction.  Because  of  the 
configuration  of  our  test  range  and  structure,  it  was  most  straightforward  to  determine  the 
absolute  accuracy  of  the  GPS  heading  indication. 


RELATIVE  PITCH 


Equipment  and  Procedure 

To  evaluate  absolute  heading  accuracy,  the  tree  structure  is  set  up  in  a  known  orientation 
as  follows.  The  structure  is  centered  on  the  base  point.  A  survey  pole  mounted  in  a 
tripod  is  centered  on  one  of  the  control  points  using  a  plumb  line.  The  structure  is  then 
rotated  until  the  laser  image  appears  on  the  survey  pole.  GPS  attitude  solution  data  is 
collected  for  approximately  15  minutes.  The  procedure  is  repeated  for  other  control 
points.  Errors  in  positioning  the  structure,  the  survey  poles,  and  aligning  the  laser,  are 
estimated  to  contribute  no  more  than  3  cm  of  systematic  error.  The  survey  error 
contributions  are  negligible  in  comparison.  Thus,  over  the  survey  baselines  of  30-70  m 
the  total  error  in  the  "known"  orientation  of  the  structure  is  less  than  0.06  degrees. 
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Since  the  laser  is  not  precisely  aligned  with  the  electrical  centerline  of  the  GPS  antenna 
array,  a  calibration  must  be  performed  each  time  it  is  mounted.  The  actual  misalignment 
involves  both  an  offset  and  a  rotation  of  the  laser  boresite  from  the  structure  centerline,  as 
shown  in  Figure  5.22.  However,  because  the  laser  image  is  viewed  from  at  least  30 
meters  away,  only  the  rotational  misalignment  is  important 


Target 


5.22.  Laser  Misalignment  from  Antenna  Centerline 


Results 

The  absolute  accuracy  of  the  receiver  was  tested  by  aligning  the  structure  with  survey 
poles  at  the  control  points  using  a  laser  as  described  in  Section  5.5.  These  tests  were 
performed  on  May  10,  12,  18,  and  20,  1993.  Data  was  collected  with  the  structure 
pointed  at  various  sites. 

The  laser  misalignment  was  determined  by  computing  the  difference  between  the  known 
heading  and  the  mean  measured  heading,  and  a  correction  based  on  this  misalignment 
was  applied  to  the  correct  the  measured  heading.  Details  of  the  laser  misalignment 
calibration  are  found  in  [Axelrad  and  Chesley,  1993].  Excluding  Point  A,  for  which  the 
absolute  heading  was  not  reliable,  the  best  agreement  was  ±0.07  degrees.  The  poorest 
agreement  for  the  tests  was  ±0.20  degrees.  Changes  in  satellite  geometry,  signal  strength, 
multipath  effects,  and  integer  bumps  account  for  the  differences  between  sessions. 
Intege^  bumps  will  be  discussed  in  the  next  subsection. 

Typical  heading  results  from  the  pointing  tests  are  shown  in  Figure  5.23  when  the 
structure  was  aligned  with  points  B  and  C.  The  mean  measured  heading  differs  from  the 
tnie  heading  by  -0.16  degrees  for  Point  B  and  by  +0.20  degrees  for  Point  C.  Figure  5.24 
summarizes  the  heading  errors  for  all  the  absolute  pointing  tests. 
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Fig.  5.23.  Adjusted  Heading  Results  for  Tree  Structure,  12  May  1993 
Heading  Error  Magnitude(deg) 
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Fig.  5.24.  Summary  of  Adjusted  Heading  Errors  for  Absolute  Tests. 
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We  compared  the  real  time  heading  reported  by  the  receiver  to  the  known  value  for  a 
nearly  orthogonal  antenna  array  with  baselines  of  1 .4  m.  The  mean  error  for  each  test 
was  0.2  deg  or  less.  The  standard  deviation  was  approximately  0.1  deg. 

The  accuracy  of  the  receiver  is  limited  primarily  by  systematic  errors,  which  are  probably 
due  to  multipath  and  antenna  phase  center  motion  as  described  above.  For  these  tests,  the 
reference  was  provided  by  a  diode  laser  mounted  to  the  structure.  This  laser  can  be 
thought  of  as  representative  of  a  spacecraft  instrument  which  is  to  be  pointed  at  a  target. 
The  offset  or  misalignment  of  the  GPS  antenna  array  centerline  to  the  boiesite  of  interest 
must  be  calibrated.  In  our  tests,  a  new  calibration  was  required  for  each  day.  However, 
this  is  very  impractical.  For  spacecraft  applications  and  for  our  future  ground  tests,  we 
plan  to  perform  this  calibration  in  conjunction  with  the  structure  baseline  self  survey. 


5.10  Degraded  Test  Structure 

An  alternate  antenna  configuration  was  used  to  compare  pointing  accuracy  results  for 
different  antenna  baseline  configurations.  If  the  baselines  are  orthogonal,  then  attitude 
accuracy  is  enhanced.  If  the  baselines  are  coincident,  then  attitude  determination  is 
degraded.  The  objective  of  this  test  was  to  gain  a  qualitative  sense  for  how  much 
degradation  could  be  expected  for  the  tree  structure. 

Equipment  and  Procedure 

The  arms  of  the  tree  structure  were  modified  so  that  all  the  slave  antennas  were  roughly 
in  the  same  direction  relative  to  the  master  and  all  four  antennas  are  more  coplanar.  This 
was  done  to  decrease  the  sensitivity  of  the  antenna  geometry  to  changes  in  heading.  So 
rather  than  being  orthogonal,  the  maximum  angular  separation  of  antenna  baselines  was 
15  degrees.  Attitude  data  were  then  collected  in  using  the  same  procedure  used  for  the 
absolute  accuracy  tests. 

Results 

Test  results  using  the  structure  with  poor  geometry  confirm  that  antenna  placement  is  an 
important  limiting  factor  in  attitude  solution  accuracy.  Not  only  should  the  baselines  be 
as  long  as  practically  possible,  but  they  should  also  be  configured  as  nearly  orthogonal  as 
possible.  Sample  test  results  for  the  geometrically  degraded  structure  are  shown  in  Figure 
5.25. 
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5.25.  Adjusted  Heading  Results  for  Degraded  Test  Structure,  2  Jun  1993 
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The  standard  deviation  in  heading  for  the  degraded  structure  is  approximately  double  that 
for  the  optimal  structure,  and  the  heading  solution  suffers  many  jumps  due  to  integer 
bumps  in  the  solution  algorithm.  This  indicates  that  the  integer  number  of  wavelengths 
between  the  antennas  was  not  well  resolved  due  to  the  poor  structure  geometry,  and  limits 
the  usefulness  of  this  structure  geometry  for  determining  accurate  attitude  results. 

The  key  to  accurate  GPS  based  attitude  determination  is  correctly  determining  the  integer 
ambiguities.  The  receiver  Trimble  automatically  monitors  the  validity  of  these  integers 
and  corrects  them  if  that  becomes  necessary.  A  faulty  integer  bump  would  temporarily 
degrade  the  attitude  solution  until  the  monitoring  algorithm  detected  it  and  returned  it  to 
the  correct  value.  An  example  of  this  can  be  seen  in  Figure  5:26.  Note  the  downward 
jump  in  the  yaw  value  at  about  3  minutes  into  the  run  corresponding  to  an  integer  bump. 
The  integer  is  bumped  back  again  after  about  IS  seconds,  and  the  heading  returns  toward 
the  correct  value.  The  entire  episode  lasts  about  30  seconds,  and  the  largest  excursion  is 
approximately  0.5  degrees.  These  effects  tend  to  average  out  for  long  static  runs,  but  they 
could  introduce  problems  in  critical  real  time  applications. 


TIME  (min) 


Heading  :  Mean  =  354. 1084  Sid  =  0. 15  deg 

ACTUAL  HEADING:  354.1864  deg 

Fig.  5.26  Adjusted  Heading  Results  for  Tree  Structure 
Showing  Incorrect  Intweger  Bump 

5.11  Dynamic  Tests 

Initial  tests  have  been  conducted  on  the  roof  of  the  Engineering  Center  using  the  one- 
armed  spinner.  The  objective  of  these  tests  was  to  gain  preliminary  experience  using  the 
receiver  in  a  dynamic  environment.  The  single  baseline  of  the  one-armed  spinner  was 
surveyed,  and  differential  phase  data  has  been  collected  while  the  platform  was  rotating. 
At  the  time  of  this  writing,  no  attitude  solutions  have  been  generated  using  the  spinning 
platform. 
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6.  Future  Work 

In  this  first  year  project,  good  progress  was  made  in  establishing  the  tools  and  expertise 
for  research  in  GPS  based  spacecraft  attitude  determination.  In  the  process,  we  have 
identified  several  key  areas  which  we  plan  to  pursue.  The  objective  of  this  continuing 
research  is  to  develop  new  algorithms  which  are  more  accurate,  reliable,  and  robust  for 
real  time,  on-orbit  operation. 

The  first  task  ahead  is  to  validate  the  existing  algorithms  using  data  collected  on-orbit. 
We  have  just  received  a  set  of  raw  data  files  collected  from  the  RADCAL  satellite.  Some 
preliminary  results  based  on  this  data  were  presented  by  C.  Cohen  and  G.  Lightsey 
[1993]  at  the  recent  ION  conference  in  Salt  Lake  City.  There  are  several  challenges  to 
working  with  the  RADCAL  data.  Unfortunately,  a  reliable  static  survey  of  the  baselines 
was  not  performed  prior  to  launch;  thus,  the  baselines  and  line  biases  must  be  determined 
from  the  on-orbit  data.  Second,  the  multipath  environment  on  the  spacecraft  is  expected 
to  be  poor  because  of  the  many  structures  located  near  the  GPS  antennas.  To  try  to 
reduce  this  effect,  the  antennas  were  canted  away  from  the  vehicle  centerline.  However, 
this  means  that  there  may  not  be  as  many  satellites  visible  to  all  antennas,  as  in  the  ideal 
case  where  the  antennas  are  oriented  opposite  nadir.  Finally,  the  geometry  is  also  not 
ideal  in  that  the  three  antenna  baselines  are  coplanar.  Of  course,  these  are  all  realistic 
considerations  in  implementation  of  GPS  for  onboard  attitude  determination;  thus  they 
provide  a  good  test  for  our  algorithms. 

Our  planned  approach  to  processing  the  RADCAL  data  is  as  follows: 

•  Process  the  pseudorange  data  through  the  navigation  filter  to  solve  for  the  vehicle  orbit. 

•  Store  the  computed  satellite-to-GPS  line  of  sight  vectors  to  a  file. 

•  Batch  process  the  phase  data  for  the  entire  5  hour  pass  to  solve  for  the  baseline 
coordinates,  line  biases,  and  integer  ambiguities. 

•  Repair  observed  cycle  slips  in  the  data. 

•  Fix  the  baselines  and  line  biases  and  pass  the  phase  data  through  the  attitude  filters  to 
solve  for  the  vehicle  orientation. 

Unfortunately,  to  our  knowledge,  there  is  no  external  attitude  reference  with  which  to 
compare  our  results.  However,  as  a  check,  we  do  plan  to  compare  both  intermediate  and 
final  results  with  those  being  generated  by  Lightsey  at  Stanford  University. 

Once  we  are  able  to  achieve  reliable  results  with  this  multi-stage  process  we  will  begin 
investigating  enhancements  to  the  filter  to  permit  real  time  estimation  of  the  baseline 
configuration  and  line  biases.  Improvements  in  our  current  method  are  needed  to  perform 
the  integer  ambiguity  resolution  and  both  detection  and  correction  of  cycle  slips.  The 
new  method  can  be  evaluated  with  both  the  RADCAL  data  and  simulated  data.  We  will 
also  continue  to  explore  other  opportunities  for  obtaining  data  collected  on-orbit  to 
further  validate  this  work. 
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The  computer  simulation  needs  to  be  upgraded  to  account  for  alternate  vehicle  dynamics 
such  as  an  inertially  pointing  spacecraft,  and  more  realistic  baseline  configurations  and 
visibilities.  We  are  also  interested  in  investigating  the  effects  of  vehicle  flexure  on  the 
performance  of  the  existing  algorithms  and  enhancing  them  to  permit  some  changes  in 
the  structure  or  antenna  configuration. 

In  the  ground  test  area  we  are  continuing  to  work  with  the  data  collected  on  the  one 
armed  spinner.  We  hope  to  extract  some  useful  information  from  these  tests  and  to  use 
the  experience  gained  in  this  experiment  to  guide  the  upgrade  of  the  test  bed.  We  are 
planning  to  improve  the  dynamic  test  platform  and  to  conduct  more  systematic  tests  in  a 
similar  configuration.  In  particular,  the  new  structure  will  be  considerably  stiffer  than 
those  used  thus  far,  and  an  external  "truth”  reference  will  be  provided  by  a  platform 
encoder. 

One  feature  that  was  evident  in  both  the  simulation  studies  and  ground  tests  is  the 
detrimental  effect  of  multipath  on  the  attitude  solutions.  We  propose  to  conduct  an 
experiment  using  narrow  beam  antennas  (similar  to  those  used  on  the  codeless  attitude 
determination  system  designed  by  MacDoran  and  Snyder)  with  the  Vector  receiver.  The 
performance  of  the  receiver  with  the  narrow  beam  versus  standard  patch  antennas  can  be 
compared  in  a  ground  test  on  Engineering  Center  roof.  Another  possibility  is  to  try  to 
accurately  simulate  the  performance  of  such  an  antenna  on  orbit  to  evaluate  the  potential 
benefits  for  spacecraft  applications. 

During  the  next  year  we  hope  to  continue  to  learn  more  about  GPS  performance  both  on 
the  ground  and  in  orbit,  and  to  develop  reliable  attitude  determination  algorithms  which 
could  be  tested  in  real-time  onboard  a  near  Earth  spacecraft.  If  we  are  successful,  the 
next  step  will  be  perform  closed-loop  onboard  attitude  control  based  on  GPS  feedback. 
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Appendices 


A.  Calculating  a  Transformation  Matrix  from  Quaternions  and  Euler  angles 


The  transformation  matrix,  C,  can  be  expressed  in  terms  of  Euler  angles.  Let  6 1  be  the 
rotation  about  the  1  -axis,  62  be  the  rotation  about  the  2-axis,  and  63  be  the  rotation  about 
the  3-axis.  Perform  these  rotations  in  3-2-1  order  where  63  is  the  first  rotation,  62  is  the 
second  rotation,  and  61  is  third  rotation.  This  gives 


C32i(0i,fl2’^)  = 


cos^cosdj  cos^2si''%  -sin02 

sindisind2COsd3  -  cosdjsin^  sin^sind2sin03 -cosdicos03  sindicos02 
[cos0isind2cosd3  +  sind|Sin&3  cos0}Sin02sin^  -  sin0]Cos03  cos0icos02 


(A-1) 


The  transformation  matrix,  C,  can  also  be  expressed  in  terms  of  quaternions  as  follows. 


C(q)  =  [ql  -  )l  +  Iqq^  +  Iq^q"" 


or. 


C(q)  = 


q^-ql-Ql-^ql 

2(^2^1-9493) 

2(9391+949:) 


2(9i92  +  9493) 
“9l  +92  “93  +94 
2(9392  “  949i) 


2(9i93  -9492) 

2(9:93 +  949i) 

2  2  2  2 
-9i  “9:  +93  +94 


'91' 

where  q  = 

92 

and  q  = 

9 

.^4. 

.^3. 

(A-2) 
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B.  Dynamics  Matrices 

The  non-zero  elements  of  the  dynamics  matrix  derived  from  the  linearized  angle 
equations  of  motion  (4.20)  and  used  in  the  Euler  angle  filter  equation  (4.54)  are 

^14  ~  ^25  ~  ^36  “  ^ 

F,,=[-nK-{h-hW]lh 

F,,=[-nK-A{I,-iy]II, 

F,,={K-{h-^h-h)n]lh 


The  elements  of  the  dynamics  matrix  derived  from  the  quaternion  equations  of  motion 
(4.19)  and  used  in  the  quaternion  filter  in  equati«.  i  (4.54)  are 


Fn=  0 

=  1(0)3+/) 

^13  = 

^14  =  1^4 

Fi5  =  -T93 
^16  =  1^2 


Fjl  =  ~2^®3  +/) 

F22=  0 

^23= 

^24  =  7^3 

^25=  1^4 

F26  =  “T9i 


F3,  -  -^0), 

F32  —  “"7  0)2 
F33  =  -■ ^(0)3  +/) 
^34  =  -i9l 
^35  =  -2  92 
^36  =  -2  93 


(B-2) 


F4,  =  -12^,fl‘(2^,^293  +9294  -9394) 
F42  =  -12/:,fl^(2<7,^294  +  9^93  -9394) 
F43  =  -l2K^Q^{-2q^q^q  +  q^q2  -qlq4) 

F44=  0 

F45  =  KjOHj  —h^  1 1\ 

^46=  ^l“2 
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/Tjl  —  — 6^2^  (^3  ~  ^92^3  ""  ^^3  ) 

F52  =  -6K.Q^{q4  -  ^qxq^q-i  - (>qiq4  - '^qlq*) 

F53  =  -eK^i^'iqx  -  ^^q^qA  -  6^i^3  “  ) 

F54  =  ^  ^  ^2 

^55=  0 
^56=  ^2“l 

F^,  =  -6K,QHq2  -  2q.ql  -  Iql) 

F^2  =  -6KjQ.^{qi  -4^2^3^4  -6^i^2  -29j93) 

F63  =  -eK-i^^i-q^-Aq^q^q-^  +%3^4  +2^294) 

F(^=  ^130)2 

^65= 

^66=  0 

where  .  K2  =  ,  A:3  =  --  . 

A  ^2  ^3 

and  0'=\i/r^. 

The  non-zero  elements  of  the  dynamics  matrix  derived  from  the  orbit  equations  of  motion 
(4. 1 )  and  used  in  the  navigation  filter  in  equation  (4.54)  are 
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The  wooden  pallets  were  used  for  tests  where  antenna  location  and  orientation  had  to  be 
maintained  precisely  from  one  day  to  the  next.  Two  wooden  pallets  were  placed 
approximately  3  m  apart  on  the  flat  roof  of  a  low  building  at  Table  Mountain.  This 
eliminated  any  structural  flexure  that  may  have  been  present  in  the  other  test  structures. 
The  master  antenna  was  mounted  to  one  of  the  pallets  with  C-clamps,  and  two  slave 
antennas  were  mounted  to  the  other  pallet.  Figure  5.7  shows  the  pallet  configuration  at 
Table  Mountain. 


Fig.  5.7  Wooden  Pallet  Structure 
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The  one-armed  spinner  structure  was  designed  to  mount  on  an  existing  spinning  platform. 
The  structure  consists  of  a  single  master-slave  antenna  pair  mounted  on  a  PVC  center 
pipe  and  a  single  1  meter  long  PVC  pipe  arm.  The  arm  was  fitted  to  a  cylinder  mounted 
on  the  spin  table.  The  table  is  driven  by  a  DC  motor  which  has  a  variable  power  supply 
that  can  be  used  to  regulate  spin  rate.  The  receiver,  power  supply,  and  data  collection 
computer  are  all  stored  inside  the  cylinder  during  testing.  The  configuration  is  shown  in 
Figure  5.6. 


Fig.  5.6  One-Armed  Spinner  Structure 
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5.5  Ground  Testing  Conflguration 

Ground  testing  of  the  Trimble  Vector  receiver  was  conducted  under  both  static  and 
dynamic  (spinning)  conditions.  The  tests  focused  on  measurement  quality,  static  self 
survey,  static  attitude  solution  accuracy,  and  dynamic  performance  on  a  spinning 
platform.  Before  addressing  the  specific  experiments,  the  test  beds  will  be  described. 
Each  test  bed  consisted  of  a  site,  an  antenna  platform,  and  known  reference  points  to 
verify  accuracy. 

Test  locations 

Two  locations  were  used  for  ground  testing  -the  roof  of  the  Engineering  Center,  and 
Table  Mountain.  The  roof  of  the  Engineering  Center  at  the  University  of  Colorado 
consists  of  a  flat  area  approximately  60  x  100  m  (see  Figure  5.4).  The  view  of  the  sky  is 
partially  obstructed  to  the  southwest  by  the  tower  of  the  Engineering  Center  and  adjacent 
buildings.  Many  antennas  and  metal  structures  are  co-located  on  the  roof,  making  it  a 
very  noisy  environment  for  multipath.  (However,  it  was  quite  convenient.) 


Fig.  5.4.  Engineering  Center  Roof  and  x  Structure 


The  other  test  location  used  was  Table  Mountain,  shown  in  Figure  5.5.  Table  Mountain 
is  a  large  flat  mesa  north  of  Boulder,  Colorado,  controlled  by  the  U.S.  Department  of 
Commerce.  The  area  used  was  largely  free  of  multipath  reflectors  and  contained 
surveyed  reference  points. 
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Fig.  5.5  Table  Mountain  Site  and  Tree  Structure 


Antenna  Platforms 

Four  antenna  platforms  were  used  for  receiver  testing.  Each  platform  was  tailored  for  a 
particular  purpose  in  the  testing  program.  The  platforms  will  be  referred  to  as  the  x 
structure,  the  tree  structure,  the  one-armed  spinner,  and  the  wooden  pallets. 

The  X  structure  is  constructed  from  two  aluminum  cross  members,  each  1-m  long  as 
shown  in  Figure  5.4.  The  master  antenna  is  raised  on  spacers  that  are  approximately 
10  cm  long  to  prevent  all  four  antennas  from  being  coplanar.  A  small  diode  laser  can  be 
mounted  at  the  center  of  the  two  cross  members  for  absolute  accuracy  tests.  The  x 
structure  was  used  in  early  tests,  but  it  suffered  from  several  shortcomings  —the  baselines 
were  too  short;  the  antennas  were  too  close  to  coplanar;  the  arms  were  too  flexible;  and 
the  structure  itself  was  a  multipath  reflector. 

Tiie  tree  structure  shown  in  Figure  5.5  was  designed  to  address  the  problems  with  the  x 
structure.  It  has  a  center  mounting  pole  for  the  master  antenna  and  three  one-meter  arms 
for  the  slave  antennas.  The  center  pole  and  baseline  2  define  the  body  coordinate  system. 
The  arms  can  be  moved  to  change  the  relative  height  and  orientation  of  the  slave 
antennas.  The  structure  is  fabricated  using  Schedule  40  PVC  pipe.  This  material  was 
chosen  because  it  was  thought  to  have  several  advantages  [Comp,  et  al.,  1993]  -  greater 
rigidity,  preventing  significant  flexure  during  the  tests;  lighter  weight,  permitting  easy 
relocation  and  storage;  and  transparency  to  L-band  GPS  signals.  The  small  diode  laser 
was  mounted  below  the  baseline  2  antenna  for  absolute  accuracy  tests.  Note  that  the 
three  antenna  baselines  shown  in  the  photograph  are  nearly  orthogonal  with  the  master 
antenna  at  the  vertex  of  the  orthogonal  set  This  nearly  orthogonal  configuration  provides 
the  best  sensitivity  to  changes  in  attitude  for  all  three  axes.  Field  experience  showed  that 
the  tree  structure  was  not  as  rigid  as  originally  hoped,  thus  suggesting  the  area  of  attitude 
determination  for  non-rigid  structures  as  an  area  for  future  research. 
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